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IV 


ABSTRACT 


This  thesis  presents  the  most  complete  climatology  of  radio  wave  (VHF  and 
higher  frequencies)  tropospheric  ducting  within  Turkey  and  the  surrounding  areas 
ever  published.  This  ducting  can  result  in  greatly  extended  propagation  ranges. 
Atmospheric  data  from  radiosonde  upper  air  observations  were  collected  for  a 
five-year  period  at  stations  in  Adana,  Ankara,  Diyarbakir,  Isparta,  Istanbul,  Izmir, 
Athens,  Bucharest,  Crete  Island,  and  Tuapse.  The  upper  level  and  surface 
ducting  yearly  average  frequencies  for  these  stations  were  59.46%,  34.22%, 
37.53%,  38.15%,  50.67%,  45.69%,  64.69%,  30.24%,  65.59%,  and  21.05%, 
respectively.  Ducts  occurred  mostly  in  summer.  These  frequencies  were  much 
higher  than  previously  published  climatologies.  However,  many  of  these  ducts 
were  too  shallow  to  be  operationally  significant.  Atmospheric  pressure  was 
correlated  with  duct  heights  for  some  stations,  but  sometimes  the  correlation  was 
positive  while  other  times  it  was  negative.  In  most  cases,  no  correlation  was 
found  between  duct  height  and  wind  speed  or  direction.  Ducts  tended  to  be 
more  common  in  coastal  regions.  Regions  away  from  the  coast  usually  have 
topography  limits  the  effects  of  the  ducts  away  from  the  coast,  and  these  inland 
regions  tended  to  have  less  ducting  occurrence  Therefore,  these  results  are  most 
applicable  to  the  coastal  regions. 
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I.  INTRODUCTION 


The  purpose  of  this  study  is  to  find  ducting  conditions  and  to  quantify 
patterns  and  trends  of  them  within  Turkey  and  surrounding  areas,  and  to 
understand  how  ducting  can  be  effective  on  the  performance  of  communication 
systems. 

Electronic  Warfare  (EW)  plays  an  important  role  in  the  overall 
scheme  of  military  strategy.  Information  superiority  plays  a  vital  role 
in  the  operation  area.  A  key  role  of  EW  is  to  be  able  to  exploit 
electromagnetic  (EM)  missions  from  a  hostile  source,  turn  that 
information  into  intelligence  on  an  adversary’s  order  of  battle, 
intentions,  and  capabilities,  and  to  render  communications  and 
weapons  systems  useless,  all  while  protecting  their  own  system. 
(Cairns-McFeeters,  1992) 

The  key  matter  is  how  the  systems  will  be  used  effectively  to  send  and 
receive  signals.  Many  things  can  affect  the  electromagnetic  radiation  in  the 
earth’s  atmosphere.  The  ducting  is  one  of  the  reasons  for  the  causes  of 
anomalous  propagation.  Through  ducting,  radio  frequency  detection  ranges  can 
be  extended.  Understanding  ducting  enables  one  to  take  advantage  of  the 
environment  in  EW. 

Anomalous  propagation  occurs  frequently  and  therefore  must  be 
factored  into  any  warfare  commander’s  decision-making  process. 
Weather  forecasters  supporting  military  operations  should 
understand  the  effects  of  atmospheric  variables  on  radio  and  radar 
performance.  They  should  collect,  analyze,  and  evaluate  available 
data  to  describe  the  existing  atmospheric  propagation  condition  and 
to  understand  how  it  will  change  over  time.  They  should  provide 
radar  operators  and  operational  commanders  with  this  information 
on  a  timely  basis.  (Ford,  2005) 

Determination  of  refractive  conditions  can  help  to  predict  the  following: 

•  How  far  away  enemy  sensors  will  be  able  to  detect  friendly  assets 
(missiles,  aircraft,  ships,  boats), 
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•  How  far  away  friendly  sensors  will  be  able  to  detect  enemy  assets, 
and  where  reconnaissance  platforms  need  to  be  placed  to  exploit 
the  atmosphere  and  extend  friendly  detection  capability, 

•  What  altitudes  provides  aircraft  with  the  greatest  weapons  stand-off 
range,  and 

•  From  what  direction  friendly  assets  should  approach  the  enemy  to 
give  the  greatest  detection  advantage  (Ford,  2005). 

This  thesis  includes  an  overview  of  EM  propagation  through  the 
atmosphere,  discussion  of  environmental  factors  that  affect  EM  propagation,  and 
finds  ducting  conditions  for  ten  different  locations.  It  concerns  radio  waves  at 
VHF  and  higher  frequencies  and  the  space  wave  mode  of  transmission. 

Atmospheric  data  was  collected  from  ten  different  radiosonde  stations. 
The  stations  examined  are  located  in  Adana,  Ankara,  Diyarbakir,  Isparta, 
Istanbul,  Izmir,  Athens,  Bucharest,  Crete  Island,  and  Tuapse.  To  find  refractivity 
profiles  for  each  radiosonde,  a  MATLAB  program  was  written.  This  program 
reads  radiosondes,  calculates  modified  refractivity,  finds  the  first  three  ducts  in 
the  profile  if  they  exist,  and  finally  does  some  statistics  to  evaluate  the  ducting 
effect  on  EM  propagation. 

This  thesis  adds  power  to  operators  and  commanders’  prediction  ability. 
Knowing  how  to  use  the  ducting  climatology  can  extend  commanders  and 
operators’  line  of  sight.  Communication  is  like  breathe:  you  have  never  known  its 
importance  until  it  is  gone. 
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II.  ELECTROMAGNETIC  WAVES  AND  ATMOSPHERIC 

PROPAGATION 


A.  ELECTROMAGNETIC  WAVE 

Electrons  have  negative  charges  and  protons  have  positive  charges.  If  a 
particle  is  charged,  it  produces  an  electric  field  in  all  directions.  While  particles  of 
the  same  charge  repel,  oppositely  charged  particles  attract  each  other.  If  the 
particle  moves,  the  electric  field  changes.  The  electrical  potential  difference 
between  two  locations  is  called  voltage. 

If  a  charged  particles  moves,  it  creates  magnetic  field.  If  an  electron 
moves  through  in  a  wire,  the  magnetic  field  encircles  the  wire.  Its  direction 
depends  on  the  direction  of  electron  movement.  If  a  charged  particle  moves,  it 
produces  both  an  electric  field  (because  the  particle  is  charged)  and  a  magnetic 
field  (because  the  particle  is  moving).  An  accelerating  charged  particle  creates 
an  Electromagnetic  (EM)  wave. 

Figure  2  shows  the  creation  of  EM  waves.  The  magnetic  and  electric 
fields  of  an  electromagnetic  wave  are  perpendicular  to  each  other  and  to  the 
direction  of  the  wave.  The  motion  of  the  particle  is  changing,  the  electric  field  is 
changing  and  the  magnetic  field  is  changing.  The  changing  electric  field  creates 
a  new  magnetic  field,  and  the  changing  magnetic  field  produces  a  new  electric 
field.  The  collapsing  and  regeneration  of  the  electric  and  magnetic  fields  is  what 
allows  EM  radiation  to  propagate  (Guest,  2005). 


z 


Figure  1 .  Electromagnetic  (EM)  Wave  Generated  by  an  Accelerating  Electron 

[From  (Guest,  2005)] 
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B.  ATMOSPHERIC  PROPAGATION 


Electromagnetic  wave  propagation  is  affected  by  the  properties  of  the 
medium  located  between  the  transmitting  antenna  and  the  receiving  antenna.  In 
most  situations  radio  waves  in  free  space  travel  in  straight  lines  (an  exception 
being  gravitational  bending  as  predicted  by  General  Relativity  theory),  while  radio 
waves  traveling  within  the  earth's  atmosphere  are  affected  by  varying  conditions. 
The  influence  exerted  on  radio  waves  by  the  earth's  atmosphere  adds  many  new 
factors  to  complicate  what  at  first  seems  to  be  a  relatively  simple  problem. 
These  complications  are  due  to  a  lack  of  uniformity  within  the  earth's 
atmosphere.  Atmospheric  conditions  vary  with  changes  in  height,  geographical 
location,  and  even  with  changes  in  time  (day,  night,  season,  and  year). 
Knowledge  of  the  composition  of  the  earth's  atmosphere  is  extremely  important 
for  understanding  wave  propagation  (Integrated  Publishing,  2007). 

Electromagnetic  radiation  (including  radio  waves)  can  be  reflected, 
refracted,  and  diffracted  within  the  atmosphere.  To  understand  radio 
propagation,  it  is  useful  to  consider  different  regions  of  the  atmosphere 
surrounding  the  earth  as  shown  in  the  Figure  2. 


Ionosphere 


Exosphere  - ► 

400 

Thermosphere 

^Mesosphere 
Stratosphere 
Troposphere 


Figure  2.  Regions  of  the  atmosphere 


Just  above  the  earth's  surface  is  the  troposphere,  followed  by  the  middle 
atmosphere  (or  stratosphere),  and  finally  the  ionosphere.  The  troposphere  is 
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characterized  by  relatively  high  moisture  content  and  is  affected  by  surface 
topography.  Virtually  all  weather  phenomena  take  place  in  the  troposphere.  The 
temperature  in  this  region  generally  decreases  with  altitude,  clouds  form,  and 
there  may  be  turbulence  due  to  wind  shear  or  surface  heating.  Vertical  variations 
in  temperature,  humidity,  density,  and  pressure  in  the  troposphere  can  have  a 
substantial  effect  on  the  propagation  of  radio  waves.  The  stratosphere  has  much 
more  constant  conditions  and  there  is  little  water  vapor  present.  Therefore,  the 
stratosphere  has  relatively  little  effect  on  radio  propagation.  The  ionosphere  has 
charged  particles  that  cause  upward  propagating  radio  waves  in  the  HF  and 
lower  bands  (see  Table  1  below)  to  be  refracted  back  toward  the  surface.  This  is 
the  most  important  region  of  the  atmosphere  for  long  distance  point-to-point 
communications  in  the  HF  range  (Davidson,  2003).  However,  this  thesis  is 
concerned  with  tropospheric  propagation. 


Table  1 .  Frequency  Ranges  and  Propagation  Characteristics  for 

Different  Bands 


Bands 

Frequency  Ranges 

Propagation  Characteristic 

ELF 

0.3  -  3  kHz 

Earth/ionosphere  waveguide  modes 

VLF 

3  -  30  kHz 

Earth/ionosphere  waveguide  modes 

LF 

30  -  300  kHz 

Ground  waves  (Guided  by  ground) 

MF 

0.3-3  MHz 

Ground  waves 

HF 

3  -  30  MHz 

Sky  waves  (Ionospheric  “Reflection”) 

VHF 

30  -  300  MHz 

Space  Waves 

UHF 

0.3 -3  GHz 

Space  Waves 

SHF 

3  -  30  GHz 

Space  Waves 

EHF 

30  -  300  GHz 

Space  Waves 
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How  an  EM  wave  behaves  in  the  atmosphere  is  usually  related  to  its 
frequency.  Different  frequencies  behave  differently  under  the  same  conditions. 
Table  1  shows  frequency  ranges  and  characteristics  for  different  bands.  Space 
waves  are  usually  considered  to  have  “line-of-sight”  propagation,  but  refractive 
effects  can  be  very  important  and  increase  or  decrease  the  line-of-sight  ranges. 

1.  Ground  Waves 

EM  radiation  can  propagate  in  various  modes  of  transmission.  Ground 
waves  are  radio  waves  that  travel  near  the  surface  of  the  earth.  The  ground 
wave  is  actually  composed  of  two  separate  component  waves:  surface  waves 
and  space  waves.  A  surface  wave  travels  along  the  surface  of  the  earth.  A 
space  wave  travels  over  the  surface. 

2.  Sky  Waves 

Sky  waves  are  radio  waves  that  are  reflected  (actually  refracted)  back  to 
earth  from  the  ionosphere. 


Figure  3.  Diagram  Showing  Ground  Waves  (Surface  and  Space)  and  Sky 

Waves  [From  (Davidson,  2003)] 
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III.  ATMOSPHERIC  EFFECTS  ON  PROPAGATION 


A.  ATMOSPHERIC  EFFECTS  ON  RADAR  RANGE 


In  the  atmosphere,  electromagnetic  waves  do  not  propagate  in  straight 
paths  as  they  do  in  free  space.  They  behave  differently  depending  on  the 
frequency  and  the  atmospheric  effects.  This  thesis  concerns  radio  waves  at  VHF 
and  higher  frequencies  and  the  space  wave  mode  of  transmission.  In  this  range, 
the  most  important  atmospheric  effects  happen  in  the  troposphere,  usually  in  the 
lower  part  near  the  earth’s  surface.  Space  waves  in  the  troposphere  undergo 
attenuation,  scattering,  and  refraction. 

For  radar  transmissions,  EM  energy  must  propagate  toward  a  target  and 
then  be  reflected  back  toward  the  antenna.  In  free  space,  radar  energy 
propagates  according  to  the  simple  form  of  the  radar  equation  (Skolnik,  2001 ): 


R 


max 


P,GAea 

M2  smi 


(3.1) 


where: 

Pt  =  Transmitted  power,  W 
G  =  Antenna  gain 
Ae  =  Antenna  effective  aperture,  m2 
o  =  Radar  cross  section  of  the  target,  m2 
Smin  =  Minimum  detectable  signal,  W. 

In  the  troposphere,  this  simple  form  of  the  radar  equation  does  not 
adequately  predict  the  range  performance  of  the  actual  radars.  There  are 
several  reasons  for  failure  of  the  simple  form  of  the  radar  equation.  This  thesis 
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will  focus  on  those  situations  where  and  when  the  simple  radar  equation  and  line- 
of-sight  propagation  do  not  occur  due  to  effects  caused  by  the  earth’s  surface 
and  atmosphere. 

Propagation  effects  might  extend  the  radar  range  or  reduce  it.  It  is 
important  to  consider  the  earth’s  environment  when  predicting  radar 
performance.  Scattering  (reflection  off  particles  and  inhomogenieties),  refraction 
(bending),  ducting,  diffraction,  attenuation,  external  noise  and  backscatter  from 
land,  sea  and  weather  clutter,  and  attenuation  in  rain  and  other  hydrometeors 
modify  free  space  radar  performance  and  affect  the  EM  propagation  (Skolnik, 
2001). 

B.  REFRACTION 

Unlike  free  space,  the  atmosphere  is  not  uniform.  Variations  in  pressure, 
temperature,  and  most  importantly,  humidity,  cause  electromagnetic  waves  to  be 
bent  (refraction).  A  bending  toward  the  surface  can  extend  the  radar  horizon  and 
increases  the  coverage  of  the  radar,  while  bending  away  from  the  earth  will  result 
in  decreased  ranges.  Refraction  can  also  introduce  an  error  in  the  measurement 
of  the  target  elevation  angle;  the  apparent  position  of  a  detected  target  is 
different  from  its  actual  position. 
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Figure  4.  Tracking  Errors  for  Standard  Atmosphere  [From  (Davidson,  2003)] 

Refraction  of  the  radar  waves  in  the  atmosphere  is  due  to  the  variation  of 
the  velocity  of  propagation  with  altitude.  The  index  of  refraction  (Equation  3.2)  is 
a  measure  of  velocity  of  propagation  in  a  medium  compared  to  the  velocity  in 
free  space.  EM  waves  travel  faster  in  mediums  with  relatively  lower  indices  of 
refraction.  The  index  of  refraction,  n,  for  a  particular  medium  is  determined  by 
the  magnitude  of  the  phase  velocity  in  space,  c,  (i.e.,  the  speed  of  light)  divided 
by  the  phase  speed  of  light  in  the  medium,  v. 

n  =  dv  (3.2) 

At  the  earth’s  surface,  the  value  of  the  index  of  refraction  is  around 
1.000315.  In  order  to  simplify  the  use  of  refraction  characteristics  it  is  useful  to 
define  a  quantity  called  the  refractivity,  N,  which  is  the  index  of  refraction  minus 
one  times  106.  This  results  in  easy  to  use  numbers  (such  as  315  instead  of 
1.000315).  The  refractivity,  N,  for  radio  waves  depends  on  atmospheric 
pressure,  P  (hPa),  Temperature,  T,  (K)  and  vapor  pressure,  e,  (hPa)  according 
to  the  following  formula. 

7V  =  («-l)xl06  =77.6^-5.6^  +  3.75xl05^-  (3.3) 
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Note  that  this  formula  is  valid  for  all  radio  frequencies  in  the  VHF  or  higher 
frequency  ranges.  The  index  of  refractivity  corresponds  to  N=315  at  sea  level. 

The  vertical  gradient  of  N  determines  how  EM  energy  will  propagate 
through  the  atmosphere  in  relation  to  the  horizontal.  If  dN/dz  (N  gradient)  is 
positive,  EM  energy  is  bent  upward  (subrefraction),  and  if  it  is  negative,  EM 
energy  is  bent  downward. 

To  identify  trapping  conditions  and  determine  the  occurrence  of  ducts,  it  is 
easy  to  use  a  modified  index  of  refraction,  M,  where  z  is  the  height  above  the 
surface  in  meters  (Equation  3.4): 

M(z)  =  W(z)  +  0.157Z .  (3.4) 

Subrcfraction 


Figure  5.  dN/dz  and  Refraction  Categories  [From  (Davidson,  2003)] 

When  dM/dz=0,  the  ray  curvature  equals  the  earth's  curvature.  This  is 
another  way  of  saying  that  when  the  N-gradient  is  -157  units  per  kilometer,  the 
ray  has  the  same  curvature  as  the  earth.  M  will  increase  with  height  in  the 
standard  atmosphere.  The  importance  of  M  is  that  for  a  duct  to  exist,  dM/dZ 
must  be  negative  somewhere  in  the  profile  (Davidson,  2003).  Figure  5  shows  the 
refraction  categories  according  to  the  N  gradient. 

M  and  N  profiles  can  be  compared  numerically  as  shown  in  Table  2. 
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Table  2.  Comparison  of  N  and  M  Gradients  for  Refractive  Conditions 


Refractive 

condition 

dN/dz 

(N  units/km) 

dM/dz 

(M  units/km) 

Distance  to 

Surface  Horizon 

Sub  refraction 

0  <  N 

157  <  M 

Reduced 

Normal 

-79  <  N  <  0 

78  <  M  <157 

Standard 

Superrefraction 

-157  <  N  <  -79 

0  <  M  <  78 

Increased 

Trapping 

N  <-157 

M  <  0 

Greatly  Increased 

A  trapping  layer  is  the  region  where  dN/dz  <  -157  units/km  or  dM/dz 
<  0.  The  ray  curvature  bends  downward  relative  to  the  earth.  A 
trapping  layer  that  does  not  extend  down  to  the  surface  is 
sometimes  referred  to  as  an  elevated  layer.  A  duct  is  the  region 
associated  with  the  trapping  layer;  duct  EM  energy  is  confined  and 
channeled  between  the  top  and  bottom  of  the  duct.  The  top  of  the 
duct  will  always  be  the  top  of  the  trapping  layer.  However,  the 
bottom  of  the  duct  can,  and  often  does,  extend  below  the  bottom  of 
the  trapping  layer.  A  duct  can  have  its  lower  boundary  at  the 
surface  (a  "surface  based  duct")  or  above  the  surface  (an  "elevated 
duct").  (Davidson,  2003) 

Refraction  is  the  bending  of  radio  waves.  Unlike  reflection,  there  is  no 
sudden  direction  change  in  refraction.  The  EM  wave  bends  as  it  passes  through 
regions  with  different  indices  of  refraction.  The  index  of  refraction  of  a  medium 
changes  with  the  wave’s  phase  velocity. 

(3.5) 


As  seen  from  Equation  3.5,  the  velocity  of  propagation  within  any  medium 
is  a  function  of  the  permittivity  (dielectric  constant)  of  the  medium.  The  index  of 
refraction,  n,  is  equal  where  er  is  the  relative  dielectric  constant. 
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Figure  6.  Snell’s  Law 

According  to  the  Snell's  Law,  when  a  propagating  wave  in  medium  1  with 
refractive  index  ni  impacts  medium  2  with  index  n2  at  an  incidence  angle  0i,  the 
refraction  angle  02  is  given  by  Equation  3.6  (Davidson,  2003): 

nx  sin 0X  =  n2  sin 02 .  (3-6) 

Although  in  the  atmosphere  n  does  not  change  abruptly,  Snell’s  Law  can 
still  be  used  to  predict  propagation  angles  resulting  from  changes  in  n.  From 
Snell’s  Law,  it  is  clear  that  if  the  0i  is  large  enough,  the  critical  angle  can  be 
exceeded  and  propagation  ray  paths  are  almost  parallel  to  the  trapping  layer. 
The  largest  possible  angle  of  incidence  that  still  results  in  a  refracted  ray  is  called 
the  critical  angle;  in  this  case,  the  refracted  ray  travels  along  the  boundary 
between  the  two  media.  In  addition,  the  vertical  position  of  a  target  can  be  wrong 
if  there  is  no  correction  for  atmospheric  refraction. 

C.  DUCTING 


Normally,  in  the  troposphere,  the  air  is  colder  at  higher  elevations. 
However,  sometimes  the  temperature  increases  with  elevation.  When  this 
happens,  it  is  known  as  a  temperature  inversion.  Temperature  inversions  are 
dynamically  stable  and  often  prevent  the  mixing  of  air  by  turbulence.  This  can 

trap  moisture  at  the  lower  levels  and  cause  a  strong  decrease  in  vapor  pressure, 
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e,  across  the  inversion.  The  increase  in  temperature  and  more  importantly  (for 
most  situations)  the  associated  decrease  in  vapor  pressure  (humidity)  can  cause 
radio  waves  to  bend  toward  the  earth  more  sharply  than  the  earth’s  curvature, 
and  hence  the  radiation  is  “trapped”  near  the  surface.  This  situation  associated 
with  inversion  creates  what  are  known  as  channels  or  ducts,  where  extended 
“over-the-horizon”  propagation  ranges  can  occur.  If  any  EM  wave  reaches  and 
enters  a  duct,  this  duct  acts  as  a  waveguide  to  send  signals  long  distances. 
Knowing  the  ducting  climatology  of  a  region  is  important  for  any  Electronic 
Warfare  System  that  transmits  or  receives  radio  signals  (Davidson,  2003). 


<» 

T3 

3 


< 


Temoerature 


Figure  7.  Temperature  Inversion  [From 
http://apollo.lsc.vsc.edu/classes/met130] 


The  gradient  of  the  modified  refractivity  with  height  is  the  best  method  to 
see  the  ducts.  Equation  3.3  was  about  refractivity,  N.  To  identify  the  trapping 
conditions  modified  refractivity  is  used: 


M(z)  =  W(z)  +  0. 1 57Z  =  77.6  -  5.6  +  3.75x10 

T(z)  T(z ) 
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P  =  Atmospheric  Pressure  in  milibar 
T  =  Temperature  in  Kelvin 
e  =  Vapor  Pressure  in  milibar 
Z  =  Height. 


Figure  8.  M  Relationship  for  Elevated  Trapping  Layer  and  Duct  [From 

(Davidson,  2003)] 

Figure  8  shows  different  types  of  ducts.  The  top  of  the  trapping  layer  is 
the  top  of  the  duct.  The  bottom  of  the  trapping  layer  can  be  determined  by 
drawing  a  vertical  line  from  the  top  of  the  duct  down  until  it  crosses  the  M  profile. 
If  it  is  a  surface-based  duct,  the  bottom  of  the  trapping  layer  is  the  ground. 

1.  Evaporation  Ducts 

Evaporation  ducts  most  commonly  occur  over  the  ocean  and  other  large 
bodies  of  water.  They  are  found  very  close  to  the  surface,  in  the  lowest  two 
meters  or  so.  The  air  right  at  the  surface  of  a  body  of  water  is  usually  nearly 
saturated  with  water  vapor,  which  means  that  the  relative  humidity,  RH,  is  100 
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percent.  Over  oceans,  salt  causes  the  surface  RH  to  be  98%.  Although  the 
relative  humidity  is  high  over  the  sea  where  the  air  touches,  it  usually  decreases 
rapidly  in  the  several  meters  just  above  the  sea  if  the  upper  air  is  not  saturated 
(i.e.,  foggy).  Because  of  this  rapid  decrease  in  humidity,  refractivity  also 
decreases  rapidly,  enough  to  refract  radio  waves  back  toward  the  surface.  Thus, 
a  low-lying  duct  traps  the  EM  waves.  Signals  can  propagate  along  the  surface  of 
the  water  in  an  evaporative  duct.  Evaporation  ducts  can  cause  the  radar  ranges 
for  targets  at  or  near  the  sea  surface  to  be  considerably  greater  than  the  free 
space  range;  although  this  effect  may  be  counteracted  by  the  increase  in  surface 
clutter  associated  with  stronger  evaporations  ducts  (Skolnik,  2001). 


Figure  9.  Mirage  Over  the  Water  Associated  with  an  Evparation  Duct  [From 

(Martin,  2007)] 

Figure  9  shows  an  “inferior”  mirage  caused  by  an  ocean  surface  that  is 
much  warmer  than  the  air.  Because  warmer  water  has  a  relatively  higher  vapor 
pressure,  this  situation  is  also  associated  with  strong  evaporations  ducts.  These 
mirages  associated  with  evaporation  ducts  are  best  seen  during  the  daytime, 
though  they  also  occur  at  night.  Evaporative  ducts  can  act  as  an  RF  mirror  (and 
sometimes  an  optical  mirror)  and  reflect  signals  from  the  top  of  the  duct  (Martin, 
2007). 
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2. 


Surface-based  Ducts 


The  base  of  a  surface-based  duct  is  at  earth’s  surface.  There  are  three 
types  of  surface-based  ducts.  One  of  them  is  the  evaporation  duct  discussed 
previously.  The  second  one  is  created  from  an  elevated  trapping  layer,  and  the 
third  one  is  created  from  a  surface-based  trapping  layer. 
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Figure  10.  Surface  Duct  for  (a)  Surface-based  Trapping  Layer  and  (b) 
Elevated  Trapping  Layer  [From  (Davidson,  2003)] 


As  seen  in  Figure  10,  the  value  of  M  at  the  duct  top  is  less  than  the  value 
of  M  at  the  surface. 


Figure  1 1 .  Trapped  Moisture  Over  Land  Due  to  Radiational  Cooling  and 
Inversion  Formation  Near  the  Surface  can  Create  Surface-based  Ducts 
Over  Land  [From  (Martin,  2007)] 
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3. 


Elevated  Ducts 


The  base  of  an  elevated  duct  lies  above  the  surface  of  the  earth.  This 
type  of  duct  is  similar  to  the  surface-based  duct.  The  only  difference  is  that  the 
bottom  of  duct  is  higher  than  ground  level.  This  situation  is  generally  caused  by 
temperature  inversions  or  decreases  in  humidity  (vapor  pressure)  with  elevation. 


Figure  12.  Elevated  Inversion  with  a  Base  of  620  m  is  Formed  by  Subsidence 
in  a  High-pressure  Region  and  is  Often  Associated  with  Elevated  Ducts 

[From  (Martin,  2007)] 
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IV.  SOUNDING  DATA  SETS  AND  MATLAB  CODES 


A.  SOUNDING  DATA  SETS 

The  goal  of  this  thesis  is  to  quantify  ducting  climatology  within  Turkey  and 
its  surrounding  areas.  All  the  atmospheric  upper  air  data  sets  were  acquired 
from  the  Department  of  Atmospheric  Science  at  the  University  of  Wyoming.  Five 
years  worth  of  data  were  collected  to  calculate  refractivity  for  ten  different 
stations.  Table  3  shows  these  stations  and  their  locations. 


Table  3.  Stations  and  their  Locations 


Station  Code 

Station  Identifier 

Location 

17351 

— 

Adana/Turkey 

17130 

— 

Ankara/Turkey 

17281 

— 

Diyarbakir/Turkey 

17240 

— 

Isparta/Turkey 

17062 

— 

Istanbul/Turkey 

17220 

— 

Izmir/Turkey 

16716 

LGAT 

Athens/Greece 

15420 

LRBS 

Bucharest/Romania 

16754 

LGIR 

Crete  Island/Greece 

37018 

— 

Tuapse/Russia 
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Figure  13.  Geographic  Locations  of  the  Stations 


The  first  few  lines  of  sample  data  are  shown  in  Figure  14.  It  is  a  detailed 
data  set  that  includes  all  necessary  parameters  to  calculate  modified  refractivity. 


17062  Observations  at  002  01  Jan  2005 


FEES 

hPa 

HGHT 

m 

TEMP 

C 

DWPT 

C 

EELH 

1 

MIXR 

g/kg 

VRC2 

deg 

SKNT 

knot 

THTA 

K 

THTE 

K 

THTV 

K 

1016.0 

39 

8.2 

7.1 

93 

6.26 

40 

12 

280.1 

297.4 

281.1 

1015.0 

47 

8.2 

7.1 

93 

6.26 

80 

7 

280.1 

297.4 

281.2 

1000.0 

167 

7.8 

6.8 

93 

6.23 

80 

13 

280.9 

298.2 

282.0 

Figure  14.  Sample  Data 
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B.  MATLAB  CODES 


The  MATLAB  codes  used  in  this  thesis  for  determining  duct  statistics  are 
based  on  Stahlhut’s  (2006)  routines  that  were  written  to  calculate  refractivity  in 
the  Arctic  regions.  The  MATLAB  routines  were  revised  to  input  data  from  the 
University  of  Wyoming  data  sets.  The  MATLAB  code  is  shown  in  the  Appendix. 

First,  all  soundings  are  read  by  MATLAB.  After  the  modified  refractivity  is 
calculated  for  every  height,  the  M  gradient  is  calculated.  According  to  the  M 
gradients,  the  ducting  conditions  are  determined  for  each  profile.  The  lowest 
three  ducts  are  examined  if  they  exist.  According  to  the  duct  bottom,  surface 
duct  existence  is  checked.  Then  key  duct  parameters  are  recorded  for  each 
station.  Finally,  an  Excel  sheet  that  includes  the  duct  heights,  duct  thickness, 
duct  gradient  and  other  information  related  to  the  station  is  created.  This  cycle  is 
done  for  every  sounding  for  a  year.  All  results  are  collected  then  saved  as  an 
Excel  (.xls)  file.  This  Excel  file  is  loaded  into  MATLAB.  To  find  statistics  for  one 
station  the  “statistics. m”  program  (see  Appendix,  Section  D)  is  run.  To  find 
statistics  for  a  month  the  ductformonth.m  file  is  run  before  “statistics. m.” 

To  calculate  a  modified  refractivity  at  a  given  height,  Equation  3.7  is  used. 
Equation  4.1,  gives  the  vapor  pressure  at  a  given  height  where  the  temperature, 
T,  is  in  Celsius: 

e(z)  =  6. 1 1 2e(1 7-67*r(z))/(243 • 1 2+T(z))  {4A) 


The  duct  top  and  duct  middle  are  the  local  minimum  and  local  maximum  of 
M,  respectively,  found  by  a  switch  in  the  sign  of  dM/dz.  The  duct  bottom  has  the 
same  value  of  M  as  at  the  duct  top.  It  can  be  calculated  by  dropping  a  straight 
vertical  line  from  duct  top  to  the  duct  bottom.  It  can  be  found  by  applying  simple 
geometry.  Equation  4.2  gives  the  height  of  duct  bottom  (Stahlhut,  2006). 


* 

z 


=  Zo 


M2 -m* 

M2-Mx 


(4.2) 
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V.  RESULTS 


A.  DUCT  FREQUENCY  STATISTICS 

The  objective  of  this  chapter  is  to  evaluate  the  results  gathered  from 
MATLAB.  In  this  study  it  is  important  to  know  and  understand  how  ducting 
conditions  can  effect  the  performance  of  communication  systems. 

A  transmitter  or  receiver  placed  in  a  surface  duct  will  experience 
significant  signal  enhancement,  and  long  detection  range.  Omni¬ 
directional  sources  will  experience  cylindrical  spreading  loses,  and 
narrow  radar  beams  will  view  the  duct  as  a  waveguide.  Existence 
of  trapping  layer  is  determined  by  a  temperature  increase  and/or 
humidity  decrease  with  height.  (Cairns-McFeeters,  1992) 

Table  4  shows  frequency  of  occurrence  of  ducts  during  the  five-year 
period  from  2005-2009  for  each  station.  Ductl  represents  the  lowest  ducts  that 
also  include  the  surface-based  ducts  in  Table  4.  Duct2  is  the  first  elevated  ducts 
that  are  higher  than  the  lowest  ducts.  Duct3  is  the  third  highest  ducts  in  M 
profiles. 

The  significant  result  is  that  the  probability  of  occurrence  of  a  duct  is 
relatively  high  for  all  locations.  If  we  compare  day  occurrence  with  night 
occurrence,  it  is  clear  that  there  is  not  much  differences  between  them.  For 
example,  Istanbul’s  day  occurrence  frequency  is  50.93  and  its  night  occurrence 
frequency  is  50.44.  Although  this  conclusion  is  also  true  for  other  trapping 
layers,  some  stations  have  significant  occurrence  differences  between  night  and 
day  occurrences  of  surface  ducts.  The  surface-based  duct  occurrence  in  Ankara, 
Istanbul,  and  Izmir  is  higher  during  day  hours  than  night  hours. 
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Table  4. 


Frequency  of  Ducts 


Present  occurrence  of  ducts  during  24 
hours 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

17351  Adana 

59.46 

23.38 

6.40 

16.05 

17131  Ankara 

34.22 

5.96 

0.56 

18.77 

17281  Diyarbakir 

37.53 

6.18 

0.56 

21.92 

17240  Isparta 

38.15 

8.29 

0.84 

22.00 

17062  Istanbul 

50.67 

16.54 

4.06 

17.10 

17220  Izmir 

45.59 

11.07 

1.58 

20.16 

16716  Athens 

64.69 

18.94 

3.27 

44.20 

15420  Bucharest 

30.24 

5.74 

0.83 

11.26 

16754  Crete  Island 

65.59 

23.81 

5.57 

29.45 

37018  Tuapse 

21.05 

2.33 

0.14 

13.36 

Present  occurrence  of  ducts  during  day 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

17351  Adana 

51.56 

19.34 

4.00 

23.10 

17131  Ankara 

43.53 

9.95 

0.85 

32.00 

17281  Diyarbakir 

34.79 

5.93 

0.81 

22.25 

17240  Isparta 

38.53 

9.73 

1.30 

24.72 

17062  Istanbul 

50.93 

17.45 

4.69 

25.92 

17220  Izmir 

52.85 

14.94 

2.34 

33.35 

16716  Athens 

64.22 

18.75 

2.84 

46.89 

15420  Bucharest 

28.50 

5.13 

0.61 

12.27 

16754  Crete  Island 

67.24 

26.55 

6.40 

35.89 

37018  Tuapse 

22.89 

2.53 

0.07 

15.24 

Present  occurrence  of  ducts  during  night 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

17351  Adana 

67.15 

27.32 

8.76 

9.11 

17131  Ankara 

25.10 

2.05 

0.28 

5.82 

17281  Diyarbakir 

40.18 

6.41 

0.33 

21.60 

17240  Isparta 

37.80 

6.90 

0.40 

19.41 

17062  Istanbul 

50.44 

15.65 

3.44 

8.44 

17220  Izmir 

38.46 

7.30 

0.84 

7.30 

16716  Athens 

65.23 

19.19 

3.61 

42.30 

15420  Bucharest 

31.97 

6.35 

1.05 

10.27 

16754  Crete  Island 

63.91 

20.87 

4.75 

22.66 

37018  Tuapse 

19.21 

2.13 

0.21 

11.48 
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Table  5.  Present  Occurrence  of  First  Duct  by  Month  (Surface-based  Duct  is  Included) 


Jan 

Feb 

Mar 

Apr 

May 

Jun 

Jly 

Aug 

Sep 

Oct 

Nov 

Dec 

17351  Adana 

33.44 

29.63 

43.73 

58.87 

68.92 

82.87 

91.52 

94.85 

71.72 

57.89 

42.81 

39.80 

17131  Ankara 

23.53 

25.00 

20.79 

26.44 

43.52 

44.48 

51.97 

51.80 

34.29 

30.00 

28.28 

32.89 

17281  Diyarbakir 

18.72 

27.19 

27.05 

37.45 

56.77 

53.95 

62.81 

55.14 

47.12 

24.57 

23.89 

21.90 

17240  Isparta 

20.20 

16.86 

23.31 

30.98 

45.85 

57.39 

66.54 

67.93 

48.78 

37.46 

25.17 

21.43 

17062  Istanbul 

25.00 

20.73 

29.47 

47.80 

60.98 

72.70 

81.46 

85.99 

59.93 

48.22 

40.94 

33.22 

17220  Izmir 

28.95 

29.71 

32.23 

42.37 

47.04 

61.59 

72.24 

65.25 

56.47 

38.06 

41.75 

31.40 

16716  Athens 

46.49 

44.08 

56.43 

69.51 

72.28 

83.12 

83.06 

83.77 

75.83 

63.30 

57.29 

53.21 

15420  Bucharest 

9.06 

6.16 

14.62 

26.76 

35.20 

51.38 

56.39 

55.19 

37.04 

31.72 

23.08 

14.94 

16754  Crete  Island 

44.68 

41.87 

54.31 

62.45 

69.62 

88.70 

84.85 

85.09 

75.90 

66.67 

66.67 

56.79 

37018  Tuapse 

3.42 

8.63 

8.26 

15.86 

27.65 

35.05 

47.37 

40.57 

30.47 

22.69 

12.89 

7.07 
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Figure  15.  Frequencies  of  First  Duct  by  Month  (Turkey) 
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Figure  16.  Frequencies  of  First  Duct  by  Month  (Athens,  Bucharest,  Crete  Island,  and  Tuapse) 
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Table  5,  Figure  15  and  Figure  16  all  show  the  distribution  of  first  ducts 
(surface-based  ducts  are  included)  by  month.  It  is  clear  that  the  ducts  are  most 
common  during  summer  time.  In  Adana,  there  is  greater  than  90%  probability  of 
ducts  occurring  in  July  and  August.  The  lowest  occurrence  of  first  ducts  within 
Turkey  was  in  February,  at  Isparta,  with  17  percent  of  duct  occurrence.  This 
percentage  rate  of  ducting  conditions  is  significant,  indicating  that  ducting  may  be 
a  factor  in  all  regions  of  Turkey  during  all  times  of  the  year.  The  highest 
occurrence  of  the  stations  outside  of  Turkey  were  Athens  and  Crete  Island, 
where  they  have  almost  the  same  duct  frequency  trend  because  they  are 
affected  by  the  same  climate  conditions. 

Comparing  these  results  to  AREPS  (Advanced  Refractive  Effects 
Prediction  System)  climatology  clearly  shows  that  the  percentage  of  duct 
occurrences  for  all  stations  is  higher  than  AREPS’s.  Figure  17  and  Figure  18  are 
the  sample  elevated  duct  frequency  distributions  by  months  for  Ankara  and 
Isparta  from  AREPS. 


Ankara  /Central,  Turkey 


I 1200Z 
I OOOOZ 


4  5 

Elevated  duct  -  %  occur  day  and  night 


Figure  17.  Elevated  Duct  Percent  Occurance  for  Day  and  Night  in  Ankara 
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Isparta,  Turkey 


Figure  18.  Elevated  Duct  Percent  Occurance  for  Day  and  Night  in  Isparta 

B.  DUCT  THICKNESS  AND  HEIGHT  STATISTICS 

The  following  formula,  based  on  experimental  results,  is  used  to  determine 
the  minimum  frequency  that  will  be  trapped  by  a  duct  of  given  thickness  (d,  in 
meters): 


f 


f.  = 

J  min 


A 


3.6x10  Hz 


m 


-3 

d T 


(5.1) 


Figure  19  represents  the  equations  graphically,  while  Table  6  summarizes 
duct  thickness  ranges.  Duct  thickness  statistics  of  the  stations  can  be  evaluated 
according  to  Table  6.  It  shows  the  common  frequency  ranges  and  the  required 
duct  depths  to  trap  a  signal  in  these  frequency  ranges  (Davidson,  2003). 


28 


Table  6. 


Minimum  Duct  Thickness  Required  for  Trapping  Determined  by 
Radar  Frequency  Band  [From  (Cairns-McFeeters,  1992)] 


BAND 

Frequency  Range 

Duct  Thickness  (m) 

VHF 

30  MHz -300  MHz 

122-610 

UHF 

300  MHz  -  1  GHz 

52-122 

L 

1  GHz -2  GHz 

32-52 

S 

2  GHz -4  GHz 

20-32 

C 

4  GHz -8  GHz 

12-20 

X 

8  GHz  -  12  GHz 

9-12 

K 

12  GHz -40GHz 

4-9 

0  200  60  400  120  600180  800  240  1000  300 

_ Duct  Thickness  (ft in) _ 

Figure  19.  Minimum  Trapping  Frequency  [From  (Davidson,  2003)] 

Table  7  displays  the  height  statistics  of  ducts  and  Table  8  shows  the 
thickness  of  the  first  three  and  surface  ducts.  Heights  are  above  sea  level. 
Surface  ducts  are  very  shallow  for  all  stations.  Note  from  Table  6  that  in  many 
cases,  only  the  highest  frequency  transmissions  will  be  trapped.  In  addition, 
topographic  features  such  as  hills  are  usually  higher  than  the  duct,  which  also 
would  prevent  the  normally  extended  ranges  associated  with  ducting  situations. 
Although  over  water  and  other  flat  locations  microwaves  may  have  extended 
ranges  due  to  the  ducts,  in  most  cases  these  shallow  surface  ducts  are  not 
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operationally  significant.  We  can  compare  duct  thicknesses  of  our  locations  with 
Table  6  above.  There  is  a  reverse  relationship  between  frequency  and  duct 
thickness.  As  frequency  increases,  the  required  duct  thickness  to  send  a  signal 
decreases.  The  mean  height  above  sea  level  of  first  ducts  in  Turkey  is  around 
850±100  meters  for  all  stations.  Bucharest  has  the  highest  mean  thickness  of 
surface  ducts  at  only  74.5  meters.  According  to  Table  6,  this  would  trap  radiation 
in  the  UHF  and  higher  frequency  band  ranges.  If  we  examine  the  other  layers, 
Adana  has  the  highest  mean  thickness  at  199.9  meters. 
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Table  7. 


Duct  Height  Statistics 


ADANA  (ABOVE  SEA  LEVEL  25m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

856.2 

1863.6 

2543.6 

59.1 

Lower  Quartile  (m) 

86.0 

1087.0 

1550.5 

36.0 

Median  (m) 

565.0 

1535.5 

2222.5 

45.0 

Higher  Quartile  (m) 

1346.0 

2504.5 

3552.5 

64.0 

Std.  Dev.  (m) 

970.1 

1076.6 

1239.3 

34.7 

ANKARA  (ABOVE  SEA  LEVEL  891m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

1824.9 

3169.3 

3266.7 

907.4 

Lower  Quartile  (m) 

901.0 

2651.0 

2973.0 

900.0 

Median  (m) 

920.0 

3186.0 

3381.0 

901.0 

Higher  Quartile  (m) 

2847.0 

3730.0 

3670.0 

910.0 

Std.  Dev.  (m) 

1201.3 

851.5 

679.8 

12.8 

DIYARBAKIR  (ABOVE  SEA  LEVEL  675m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

1417.2 

3267.6 

4213.9 

695.8 

Lower  Quartile  (m) 

685.0 

2362.0 

3589.0 

684.0 

Median  (m) 

706.5 

3509.5 

3971.0 

691.0 

Higher  Quartile  (m) 

1510.0 

4218.0 

4987.5 

700.0 

Std.  Dev.  (m) 

1293.9 

1413.4 

872.4 

18.4 

ISPARTA  (ABOVE  SEA  LEVEL  997m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

1759.2 

3461.2 

3502.7 

1022.7 

Lower  Quartile  (m) 

1016.0 

2960.0 

2910.5 

1007.0 

Median  (m) 

1036.0 

3507.0 

3285.0 

1016.0 

Higher  Quartile  (m) 

2706.0 

3949.0 

4030.0 

1027.0 

Std.  Dev.  (m) 

1122.6 

886.8 

753.7 

19.9 

ISTANBUL  (ABOVE  SEA  LEVEL  39m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

787.9 

1668.0 

2447.5 

63.5 

Lower  Quartile  (m) 

65.0 

941.5 

1563.0 

48.0 

Median  (m) 

419.0 

1485.0 

2373.0 

56.0 

Higher  Quartile  (m) 

1260.0 

2264.0 

3257.5 

68.0 

Std.  Dev.  (m) 

927.7 

1003.4 

1094.7 

24.4 

IZMIR  (ABOVE  SEA  LEVEL  29m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

952.9 

2062.5 

2528.6 

51.0 

Lower  Quartile  (m) 

49.0 

1292.0 

1531.0 

38.0 

Median  (m) 

174.0 

1971.0 

2470.0 

46.0 

Higher  Quartile  (m) 

1718.0 

2718.0 

3404.5 

57.0 

Std.  Dev.  (m) 

1193.0 

1056.9 

1217.1 

15.5 
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ATHENS  (ABOVE  SEA  LEVEL  15m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

505.1 

1864.6 

2747.9 

53.8 

Lower  Quartile  (m) 

39.0 

743.0 

1928.5 

32.0 

Median  (m) 

59.5 

1849.0 

2773.0 

47.0 

Higher  Quartile  (m) 

186.0 

2774.0 

3390.0 

65.0 

Std.  Dev.  (m) 

979.8 

1261.9 

1300.7 

31.3 

BUCHAREST  (ABOVE  SEA  LEVEL  91m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

953.0 

2540.9 

3443.2 

157.6 

Lower  Quartile  (m) 

154.0 

1788.0 

2858.0 

139.0 

Median  (m) 

253.5 

2474.0 

3407.0 

150.0 

Higher  Quartile  (m) 

1585.0 

3161.0 

3953.0 

167.0 

Std.  Dev.  (m) 

1153.9 

1155.8 

899.7 

33.2 

CRETE  ISLAND  (ABOVE  SEA  LEVEL  39m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

761.5 

1515.4 

2428.7 

80.1 

Lower  Quartile  (m) 

73.0 

488.0 

1270.0 

56.0 

Median  (m) 

189.5 

1193.0 

2126.0 

72.0 

Higher  Quartile  (m) 

1238.0 

2142.5 

3311.5 

95.0 

Std.  Dev.  (m) 

1021.8 

1251.7 

1513.9 

39.0 

TUAPSE  (ABOVE  SEA  LEVEL  95m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Height  (m) 

589.2 

1652.5 

1539.5 

133.7 

Lower  Quartile  (m) 

116.0 

788.0 

1034.0 

109.5 

Median  (m) 

139.0 

1566.0 

1575.0 

121.0 

Higher  Quartile  (m) 

782.0 

2207.5 

2045.0 

136.0 

Std.  Dev.  (m) 

1087.9 

990.8 

608.1 

50.5 
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Table  8. 


Duct  Thickness  Statistics 


ADANA  (ABOVE  SEA  LEVEL  25m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

147.5 

199.9 

156.6 

58.2 

Lower  Quartile  (m) 

35.0 

91.0 

78.6 

14.0 

Median  (m) 

104.3 

155.5 

119.3 

26.0 

Higher  Quartile  (m) 

212.2 

277.6 

184.9 

139.0 

Std.  Dev.  (m) 

140.4 

144.8 

121.2 

74.9 

ANKARA  (ABOVE  SEA  LEVEL  891m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

52.2 

95.1 

69.7 

16.4 

Lower  Quartile  (m) 

10.0 

59.0 

41.1 

9.0 

Median  (m) 

28.0 

87.9 

62.2 

10.0 

Higher  Quartile  (m) 

87.3 

124.7 

88.4 

58.0 

Std.  Dev.  (m) 

52.5 

49.7 

34.3 

12.8 

DIYARBAKIR  (ABOVE  SEA  LEVEL  675m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

45.2 

81.8 

79.3 

22.6 

Lower  Quartile  (m) 

10.0 

46.3 

47.8 

9.0 

Median  (m) 

26.0 

75.0 

80.6 

17.0 

Higher  Quartile  (m) 

68.8 

113.9 

112.8 

52.1 

Std.  Dev.  (m) 

45.0 

48.6 

36.2 

22.0 

ISPARTA  (ABOVE  SEA  LEVEL  997m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

51.8 

90.3 

66.3 

28.4 

Lower  Quartile  (m) 

19.0 

60.7 

39.4 

10.0 

Median  (m) 

36.2 

83.4 

50.6 

19.0 

Higher  Quartile  (m) 

75.9 

116.5 

92.5 

57.0 

Std.  Dev.  (m) 

45.7 

41.5 

40.3 

24.5 

ISTANBUL  (ABOVE  SEA  LEVEL  39m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

116.4 

153.7 

122.4 

41.8 

Lower  Quartile  (m) 

24.0 

85.9 

60.3 

9.0 

Median  (m) 

84.0 

133.5 

104.5 

18.0 

Higher  Quartile  (m) 

170.0 

191.1 

166.8 

112.6 

Std.  Dev.  (m) 

113.3 

100.2 

79.7 

69.5 

IZMIR  (ABOVE  SEA  LEVEL  29m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

70.4 

121.7 

118.0 

26.6 

Lower  Quartile  (m) 

17.0 

76.9 

53.9 

9.0 

Median  (m) 

44.0 

112.6 

90.7 

17.0 

Higher  Quartile  (m) 

105.6 

158.5 

160.7 

71.9 

Std.  Dev.  (m) 

67.3 

66.6 

80.9 

31.9 
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ATHENS  (ABOVE  SEA  LEVEL  15m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

59.3 

116.9 

99.3 

42.0 

Lower  Quartile  (m) 

24.0 

69.4 

63.0 

17.0 

Median  (m) 

40.0 

106.1 

86.9 

32.0 

Higher  Quartile  (m) 

81.6 

154.1 

123.6 

68.9 

Std.  Dev.  (m) 

54.5 

67.8 

53.8 

37.4 

BUCHAREST  (ABOVE  SEA  LEVEL  91m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

93.8 

126.5 

87.7 

74.5 

Lower  Quartile  (m) 

53.0 

79.6 

45.6 

47.2 

Median  (m) 

73.8 

111.1 

77.7 

61.0 

Higher  Quartile  (m) 

121.0 

167.0 

131.0 

109.0 

Std.  Dev.  (m) 

60.6 

67.0 

52.7 

43.5 

CRETE  ISLAND  (ABOVE  SEA  LEVEL  675m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

96.5 

138.5 

107.0 

51.8 

Lower  Quartile  (m) 

28.6 

71.0 

50.3 

19.0 

Median  (m) 

65.0 

124.8 

93.7 

34.0 

Higher  Quartile  (m) 

139.7 

188.3 

139.5 

106.0 

Std.  Dev.  (m) 

90.2 

87.6 

82.6 

59.0 

TUAPSE  (ABOVE  SEA  LEVEL  997m) 

Ductl 

Duct2 

Duct3 

Sf.  Duct 

Mean  Thickness  (m) 

76.2 

135.8 

110.9 

39.8 

Lower  Quartile  (m) 

19.0 

84.0 

55.8 

15.0 

Median  (m) 

38.0 

106.5 

105.9 

26.0 

Higher  Quartile  (m) 

96.9 

156.2 

165.9 

68.0 

Std.  Dev.  (m) 

94.0 

86.8 

75.6 

51.7 
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C.  CORRELATIONS  WITH  METEOROLOGICAL  PARAMETERS 

To  provide  guidance  on  what  factors  may  be  related  to  duct  occurrence, 
surface  pressure,  wind  speed,  and  wind  direction  were  plotted  against  duct 
height  to  determine  if  there  are  any  correlations  that  could  be  used  to  make  more 
accurate  predictions  of  duct  occurrence.  Although  these  meteorological 
parameters  are  not  directly  related  to  ducting  (which  depends  primarily  on  vertical 
variations  in  humidity),  the  hypothesis  is  that  these  parameters  may  affect  the 
humidity  profiles  which  in  turn  would  be  related  to  ducting  probabilities.  For 
example,  wind  speed  may  affect  atmospheric  mixing  of  water  vapor  (humidity) 
which  results  in  changes  in  ducting  statistics.  Wind  direction  could  be  correlated 
with  moisture  advection  (horizontal  movement),  especially  in  coastal  regions 
where  it  would  be  expected  than  air  coming  from  over  the  water  would  have 
different  humidity,  and  hence  ducting  characteristics.  Higher  surface  pressure  is 
associated  with  sinking  air  motion,  which  is  often  associated  with  the  formation  of 
an  elevated  temperature  inversion  (increase  in  temperature  with  height),  and  a 
strong  decrease  in  humidity  going  up  through  the  inversion.  It  is  well  known  that 
in  most  locations  ducting  is  correlated  with  high  pressure. 

Scatter  plots  of  these  atmospheric  parameters  and  duct  height  were 
created  for  all  stations.  Only  those  stations  with  significant  correlations  are 
displayed  in  this  section.  The  most  significant  conclusion  for  wind  speed  and 
wind  direction  is  that  they  do  not  appear  to  have  a  significant  effect  on  duct 
heights  for  most  stations,  though  some  exceptions  are  shown  below.  It  can  also 
be  seen  that,  from  plots  for  some  stations,  there  is  a  positive  relation  between 
duct  height  and  surface  pressure.  This  is  somewhat  unexpected  because  the 
sinking  air  associated  with  high  pressure  was  expected  to  push  the  duct  heights 
lower. 
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Mtop-Surface  Pressure,  ff  there  is  duel,  Adana,  Above  Sea  Level:  25m 


Figure  20.  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Adana  (if 

a  Duct  Presents) 

Figure  20  shows  that  a  higher  pressure  was  associated  with  lower  ducts  at 
Adana.  Isparta  (Figure  22)  also  has  the  same  pattern.  If  we  examine  Figure  21, 
which  is  for  surface  ducts  only,  it  is  clearly  determined  that  if  the  surface  pressure 
increases,  the  duct  height  also  increases.  This  may  be  because  the  shallower 
ducts  are  associated  with  shallow  temperature  inversions  caused  by  radiational 
cooling  (which  does  not  require  high  pressure)  while  the  deeper  ducts  may  be 
more  associated  with  the  elevated  inversions  caused  by  sinking  air. 
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Mtop-Surface  Pressure,  if  there  is  surface  duct,  Adana,  Above  Sea  Level:  25m 


Figure  21 .  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Adana  (if 

only  a  Surface  Duct  Presents) 
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Figure  22.  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Isparta  (if 

a  Duct  Presents) 
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Mtop-Wind  Speed,  if  there  is  duct,  Isparta,  Above  Sea  Level:  997m 


Figure  23.  Height  (m)  and  Wind  Speed  (knot)  Correlation  for  Isparta  (if  a  Duct 

Presents) 

For  Isparta,  there  is  a  positive  relationship  between  wind  speed  and  duct 
height  if  there  is  any  type  of  duct  (Figure  23). 

Figure  24  and  Figure  28  are  other  examples  of  a  positive  relationship 
between  surface  pressure  and  duct  height,  if  any  duct  is  present  or  if  only  a 
surface  duct  exists.  Figure  26  shows  that  if  the  wind  direction  switches  form  east 
through  north  in  Istanbul,  it  is  possible  to  see  more  numbers  of  ducts.  In 
Bucharest,  if  the  wind  speed  increases  the  duct  height  also  increases  (Figure 
27). 
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Mlop-Surface  Pressure,  If Ibere  is  duct,  Istanbul,  Above  Sea  Level  39m 


Figure  24.  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Istanbul  (if 

a  Duct  Presents) 


Mtop-Surface  Pressure,  if  there  is  surface  duct,  Istanbul.  Above  Sea  Level  39m 


Figure  25.  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Istanbul  (if 

Only  a  Surface  Duct  Presents) 
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Mtop-Wind  Direction,  if  there  is  duct,  Istanbul,  Above  Sea  Level  39m 
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Figure  26.  Height  (m)  and  Wind  Direction  (degree)  Correlation  for  Istanbul  (if  a 

Duct  Presents) 


Mtop-Wind  Speed,  if  there  is  duct,  Bucharest,  Above  Sea  Level:  91  m 


Wind  Speed 


Figure  27.  Height  (m)  and  Wind  Speed  (knot)  Correlation  for  Bucharest  (if  a 

Duct  Presents) 
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Mlop-Surface  Piessure,  If  there  is  duct.  Crete  Island,  Above  Sea  Level  39m 


Figure  28.  Height  (m)  and  Surface  Pressure  (milibar)  Correlation  for  Crete 

Island  (if  a  Duct  Presents) 
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VI.  CONCLUSION 


The  purpose  of  this  thesis  was  to  characterize  ducting  climatology  within 
Turkey  and  its  surrounding  areas.  Almost  36,500  lines  of  data,  representing 
different  elevations,  stations  and  times,  were  processed.  Refractivity  profile 
statistics  for  each  station  were  calculated.  The  general  trend  of  duct  occurrences 
and  the  relationship  between  duct  height  and  other  parameters  were  quantified. 

According  to  the  results,  the  probability  of  duct  occurrence  is  higher  in 
summer  months.  Wind  speed  and  wind  direction  do  not  have  much  effect  on 
ducting  for  most  stations.  There  is  a  positive  relationship  between  surface 
pressure  and  duct  height  at  some  stations.  This  condition  was  not  an  expected 
situation.  A  likely  explanation  is  that  very  shallow  ducts,  associated  with  surface 
radiational  cooling  or  other  local  effects,  are  often  present  during  both  high  and 
low  pressure  situations,  but  the  higher  ducts  are  associated  more  with  higher 
pressure.  Therefore,  when  considering  all  ducts,  the  average  height  is  higher  for 
high-pressure  situations.  Although  there  are  many  ducts,  they  are  not  so 
important  for  VHF  and  lower  frequencies  because  their  thickness  is  not  enough 
to  send  a  signal  in  VHF  and  lower  frequency  ranges.  These  ducts  can  be 
important  for  systems  that  use  a  microwave  frequency  range  (UHF  and  EHS). 
For  example,  radar  uses  microwave  radiation  to  detect  the  range  of  remote 
objects. 

The  plots,  tables,  and  figures  provided  in  this  thesis  can  be  a  guide  for 
Electronic  Warfare  Systems  and  radio  operators  and  their  commanders.  By  using 
knowledge  of  the  environment  and  its  affect  on  radio  transmissions,  we  can  gain 
an  advantage  over  adversaries  who  do  not  consider  this  information  in  planning 
and  carrying  out  military  operations.  It  can  also  be  helpful  for  civilian  use  of  radio 
transmissions.  The  most  significant  threat  is  losing  communication  links. 
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Knowing  the  technical  information  about  an  EW  system  and  having  the  ability  to 
predict  atmospheric  conditions  can  greatly  enhance  the  communication 
capabilities  of  Turkey  and  its  allies. 
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APPENDIX 


A.  MATLAB  PROGRAM  TO  READ  DATA 

function  loadsnd (filename) 
fid  in=fopen (filename, ' r ') ; 
if  fid_in  >  0 

sprintf (' Sounding  file  opened:  %s filename) 

end 

tic 

%  counter  for  all  data 
bigcounter  =  0; 

%This  loop  for  all  year 
while  1 

%  Skip  first  five  rows  header 
for  ihead  =  1:5; 

line  =  fgetl(fid  in); 

end 

o, 

o 

%initialize  a  line  counter 
iline=0 ; 

Q, 

O 

%This  loop  reads  each  line 
while  1 

%Read  each  line  as  a  character  string 
line=fgetl ( f id  in)  ; 

%This  loop  finds  missing  value  for  dwpt 
while  line  (22:28)  ==  ' 
line=fgetl ( f id  in); 

end 

%This  is  statement  for  reading  headers 
if  line  (1:7)  ==  'Station' 
line=fgetl ( f id  in); 
line=fgetl ( f id  in); 
station  =  str2num ( line ( 4 6 : 50 ) ) ; 
line=fgetl ( f id  in); 
date=str2num ( line ( 4  6 : 51 ) ) ; 
time=str2num (line (53:56) ) ; 
line=fgetl ( f id  in); 

%This  loop  skip  other  lines 
while  ~isempty (line)  %line(44)  ==  ':' 
line=fgetl ( f id  in); 

%Check  for  the  effective  "end-of-f ile" 
if  line  ==  -1 
break 

end 

end 

end 

%Check  for  the  end  of  one  data  set 
if  isempty (line) 
break 

end 
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%Check  for  the  effective  "end-of-f ile" 
if  line  ==  -1 
break 

end 

%Increment  line  counter.  This  is  the  index  for  each  data  array 
iline  =  iline+1; 

%Reads  data  from  file 
pres(iline)  =  str2num ( line ( 1 : 7 ) ) ; 
disp (pres (iline) ) ; 

hght(iline)  =  str2num ( line ( 8 : 14 ) ) ; 
tmpc (iline)  =  str2num ( line ( 15 : 2 1 ) ) ; 
dwpc (iline)  =  str2num ( line (22 : 2 8 ) ) ; 
relh(iline)  =  str2num (line (29 : 35) ) ; 
mixr(iline)  =  str2num ( line ( 36 : 42 ) ) ; 


if  line (48:56)  == 

drct (iline)  =  NaN; 
sknt (iline)  =  NaN; 

else 


drct (iline)  =  str2num ( line ( 47 : 4 9 ) ) ; 
sknt(iline)  =  str2num (line (54 : 56) ) ; 

end 


end  %  This  is  the  end  of  the  WHILE  Loop 


%  Calculate  theta  &  q 

[thta,q]  =  theta_q (pres, tmpc, relh) ; 

%  Calculate  M  and  N 

[M,N]  =  m  n  profile (pres, tmpc, relh, hght) ; 

%  Calculates  the  Modified  Refractivity  gradient 

Mgrad=NaN; 

for  n=l:iline-l  ; 

Mgr ad (n, 1 ) = (M (n+1 )  -  M(n))  /(hght(n+l)  -  hght(n)); 

end 

bigcounter  =  bigcounter  +  1; 

%  This  portion  finds  the  first  three  ducts,  if 

%  they  exist 

ducttop=[0  0  0] ; 

ductmiddle= [ 0  0  0]  ; 

ductbottom= [ 0  0  0] ; 

Ttop=NaN; 

Dewtop=NaN; 

Presstop=NaN; 

Mtop=NaN; 

Tmid=NaN; 

Dewmid=NaN; 

Pressmid=NaN; 

Mmid=NaN; 

Mtop2=NaN; 

Mmid2=NaN; 

Mtop3=NaN; 

Mmid3=NaN; 
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%  Find  tops  and  mids  of  ducts. . . . 

topcounter=0 ; 

midcounter=0  ; 

topindex=[0  0  0] ; 

midindex=[0  0  0]  ; 

for  m=2 : length (Mgrad) 

if  ( (Mgrad (m-1)  <  0)  &&  (Mgrad 

topcounter=topcounter+l ; 
topindex (topcounter) =m; 


end 

if  ( (Mgrad (m-1 )  >  0)  &&  (Mgrad 

midcounter=midcounter+l ; 
midindex (midcounter) =m; 


end 


(m) 


(m) 


>  0)  ) 


<  0)  ) 


end 

%  Align  these  depending  on  case . 

%  If  the  gradient  starts  negative  from  the  surface, 

%  the  subsequent  indexing  is  different  that  for  the 
%  case  if  the  gradient  starts  positive  from  the  surface, 
if  (Mgrad (1)  <  0) 

if  (topindex (1)  >  0) 

ducttop ( 1 ) =hght (topindex ( 1 ) ) ; 
ductmiddle (1) =hght (1) ; 
ductbottom ( 1 ) =hght ( 1 ) ; 

end 

if  (topindex (2)  >  0)  &&  (midindex (1)  >  0) 

ducttop (2 ) =hght (topindex (2 ) ) ; 
ductmiddle ( 2 ) =hght (midindex ( 1 ) ) ; 
if  M (midindex ( 1 ) -1 )  <  M (topindex (2 ) ) 
ductbottom (2 ) =ductmiddle (2 )  -  ... 

(M (midindex ( 1 ) ) -M (topindex (2 )))/..  . 

(M (midindex ( 1 ) ) -M (midindex ( 1 ) -1 )  )  *  .  .  . 

(hght (midindex ( 1 ) )  -  hght (midindex (1) -1) ) ; 
elseif  M (midindex ( 1 ) -1 )  >  M (topindex (2 )) &&  M (midindex ( 1 ) 

2 )  >M (midindex ( 1 ) -1 ) 

ductbottom (2 ) =hght (midindex ( 1 ) — 1 ) ; 
elseif  M (midindex ( 1 ) -2 )  >  M (topindex (2 )) &&  M (midindex ( 1 ) 

3 )  >M (midindex ( 1 ) -2 ) 

ductbottom (2 ) =hght (midindex ( 1 ) —  2 ) ; 
elseif  M (midindex ( 1 ) -1 )  >  M (topindex (2 ) ) 

counterindex  =  1; 

while  M (midindex ( 1 ) -counter index) >M (topindex (2 ) ) && . . 
(midindex (1) -counterindex) ~=1 
counterindex  =  counterindex+1 ; 


end 

ductbottom (2 ) =hght (midindex ( 1 ) -counterindex+2 )  -  ... 

(M (midindex ( 1 ) -counterindex+2 ) -M (topindex (2 ) ) ) / . 
(M (midindex ( 1 ) -counterindex+2 ) -M (midindex ( 1 ) — 

counterindex) ) * . . . 

(hght (midindex (1) -counterindex+2)  - 
hght (midindex ( 1 ) -counterindex) )  ; 
end 

Mtop2=M (topindex (2 ) ) ; 

Mmid2=M (midindex ( 1 ) ) ; 

end 
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if  (topindex(3)  >  0)  &&  (midindex (2 )  >  0) 

ducttop (3) =hght (topindex (3) ) ; 
ductmiddle ( 3 ) =hght (midindex ( 2 ) )  ; 
if  M (midindex (2 ) -1 )  <  M (topindex (3) ) 
ductbottom ( 3 ) =ductmiddle ( 3 )  -  ... 

(M (midindex (2 ) ) -M (topindex ( 3 )))/.. . 

(M (midindex (2 ) ) -M (midindex (2 ) -1 ) ) * . . . 

(hght (midindex (2 ) )  -  hght (midindex (2 ) -1 )) ; 

elseif  M (midindex (2 ) -1 )  >  M (topindex (3) ) 
counterindex  =  1; 

while  M (midindex (2 ) -counter index) >M (topindex ( 3 )  ) 
counterindex  =  counterindex+1 ; 

end 

ductbottom ( 3 ) =hght (midindex (2 ) -counterindex+2  )  -  .. 

(M (midindex (2 ) -counterindex+2 ) -M (topindex ( 3 ) ) ) / 
(M (midindex (2 ) -counterindex+2 ) -M (midindex (2 ) - 

counterindex) ) * . . . 

(hght (midindex (2) -counterindex+2)  - 
hght (midindex (2 ) -counterindex) )  ; 
end 

Mtop3=M (topindex ( 3 ) ) ; 

Mmid3=M (midindex (2 ) ) ; 

end 

Ttop=tmpc (topindex ( 1 ) ) ; 

Dewtop=dwpc (topindex (1) ) ; 

Presstop=pres (topindex ( 1 ) ) ; 

Mtop=M (topindex ( 1 ) ) ; 

Tmid=tmpc ( 1 ) ; 

Dewmid=dwpc  ( 1 ) ; 

Pressmid=pres  (1) ; 

Mmid=M ( 1 ) ; 

end 

%  These  are  for  the  cases  when  the  gradient  starts 
%  positive,  and  there  is  a  duct  present  somewhere  in 
%  the  profile 

if  (Mgrad(l)  >  0)  &  (topindex(l)  >  0) 

if  (topindex (1)  >  0)  &  (midindex (1)  >0)  &  ... 

(M (topindex ( 1 ) )  <  M(l)) 
ducttop ( 1 ) =hght (topindex ( 1 ) ) ; 
ductmiddle ( 1 ) =hght (midindex ( 1 ) ) ; 
ductbottom ( 1 ) =hght ( 1 ) ; 

end 

if  (topindex (1)  >  0)  &  (midindex (1)  >0)  &  ... 

(M (topindex ( 1 ) )  >  M(l)) 

ducttop ( 1 ) =hght (topindex ( 1 ) ) ; 
ductmiddle ( 1 ) =hght (midindex ( 1 ) ) ; 
if  M (midindex ( 1 ) -1 ) <M (topindex ( 1 ) ) 

ductbottom ( 1 ) =ductmiddle ( 1 )  -  ... 

(M (midindex ( 1 ) ) -M (topindex ( 1 )))/.. . 

(M (midindex ( 1 ) ) -M (midindex ( 1 ) -1 ) ) * . . . 

(hght (midindex ( 1 ) )  -  hght (midindex (1) -1) ) ; 
elseif  M (midindex ( 1 ) -1 ) >M (topindex ( 1 ) ) 
counterindex  =  1; 
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while  M (midindex ( 1 ) -counterindex) >M (topindex ( 1 ) ) 
counterindex  =  counterindex+1; 

end 

ductbottom ( 1 ) =hght (midindex ( 1 ) -counterindex+2  )  -  .. 

(M (midindex ( 1 ) -counterindex+2 ) -M (topindex ( 1 ) ) ) / 
(M (midindex ( 1 ) -counterindex+2 ) -M (midindex ( 1 ) — 

counterindex) ) * . . . 

( (hght (midindex ( 1 ) -counterindex+2 ) )  - 

hght (midindex ( 1 ) -counterindex) )  ; 
end 

end 

if  (topindex (2)  >  0)  &  (midindex (2)  >  0) 

ducttop (2 ) =hght (topindex (2 ) ) ; 
ductmiddle ( 2 ) =hght (midindex ( 2 ) ) ; 
ductbottom (2 ) =ductmiddle (2 )  -  ... 

(M (midindex (2 ) ) -M (topindex (2 )))/..  . 

(M (midindex (2 ) ) -M (midindex (2 ) -1 ) ) * . . . 

(hght (midindex (2 ) )  -  hght (midindex (2 ) -1 )) ; 

if  M (midindex (2 ) -1 ) <M (topindex (2  )  ) 

ductbottom (2 ) =ductmiddle (2 )  -  ... 

(M (midindex (2 ) ) -M (topindex (2 )))/.. . 

(M (midindex (2 ) ) -M (midindex (2 ) -1 ) ) * . . . 

(hght (midindex (2 ) )  -  hght (midindex (2 ) -1 )) ; 

el seif  M (midindex (2 ) -1 ) >M (topindex (2  )  ) 
counterindex  =  1; 

while  M (midindex (2 ) -counter index) >M (topindex (2 ) ) && . 
(midindex (2) -counterindex) ~=1 
counterindex  =  counterindex+1; 

end 

ductbottom (2 ) =hght (midindex (2 ) -counterindex+2 )  -  .. 

(M (midindex (2 ) -counterindex+2 ) -M (topindex (2 ) ) ) / 
(M (midindex (2 ) -counterindex+2 ) -M (midindex (2 ) - 

counterindex) ) * . . . 

( (hght (midindex (2) -counterindex+2) )  - 

hght (midindex (2 ) -counterindex) )  ; 
end 

Mtop2=M (topindex (2 ) ) ; 

Mmid2=M (midindex (2 ) ) ; 

end 

if  (topindex (3)  >  0)  &  (midindex (3)  >  0) 
ducttop (3) =hght (topindex (3) ) ; 
ductmiddle (3) =hght (midindex (3) )  ; 
ductbottom ( 3 ) =ductmiddle ( 3 )  -  ... 

(M (midindex (3) ) -M (topindex (3) ) ) / . . . 

(M (midindex ( 3 ) ) -M (midindex ( 3 ) -1 ) ) * . . . 

(hght (midindex ( 3 ) )  -  hght (midindex (3) -1) ) ; 
if  M (midindex (3) -1) <M (topindex (3) ) 

ductbottom ( 3 ) =ductmiddle ( 3 )  -  ... 

(M (midindex (3) ) -M (topindex (3) ) ) / . . . 

(M (midindex ( 3 ) ) -M (midindex ( 3 ) -1 ) ) * . . . 

(hght (midindex ( 3 ) )  -  hght (midindex (3) -1) ) ; 
el seif  M (midindex ( 3 ) -1 ) >M (topindex ( 3 ) ) 
counterindex  =  1; 
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while  M (midindex (3) -counterindex) >M (topindex (3)  ) 
counterindex  =  counterindex+1 ; 

end 

ductbottom ( 3 ) =hght (midindex ( 3 ) -counterindex+2  )  -  .. 

(M (midindex ( 3 ) -counterindex+2 ) -M (topindex ( 3 ) ) ) / 
(M (midindex ( 3 ) -counterindex+2 ) -M (midindex ( 3 ) - 

counterindex) ) * . . . 

( (hght (midindex (3) -counterindex+2) ) - 
hght (midindex (3) -counterindex) )  ; 
end 

Mtop3=M (topindex ( 3 ) ) ; 

Mmid3=M (midindex ( 3 ) )  ; 

end 

Ttop=tmpc (topindex ( 1 ) ) ; 

Dewtop=dwpc (topindex ( 1 ) )  ; 

Presstop=pres (topindex (1) )  ; 

Mtop=M (topindex ( 1 ) )  ; 

Tmid=tmpc (midindex ( 1 ) ) ; 

Dewmid=dwpc (midindex (1)  )  ; 

Pressmid=pres (midindex ( 1 )  )  ; 

Mmid=M (midindex ( 1 )  )  ; 

if  topindex (2 ) ~=0  &&  midindex (2 ) ~=0 
Mtop2=M (topindex (2 ) ) ; 

Mmid2=M (midindex (2 ) ) ; 

else 

Mtop2=NaN; 

Mmid2=NaN; 

end 

if  topindex ( 3 ) ~=0  &&  midindex ( 3 ) ~=0 
Mtop3=M (topindex ( 3 ) ) ; 

Mmid3=M (midindex ( 3 ) ) ; 

else 

Mtop3=NaN; 

Mmid3=NaN; 

end 

end 

if  ductbottom ( 1 )  <  0 

ductbottom ( 1 ) =ductmiddle ( 1 ) - [ducttop ( 1 ) -ductmiddle ( 1 ) ] ; 

end 

if  ductbottom (2 )  <  0 

ductbottom (2 ) =ductmiddle (2 ) - [ducttop (2 ) -ductmiddle (2 ) ] ; 

end 

if  ductbottom ( 1 )  <  0 
ductbottom ( 1 )  =0; 

end 

if  ductbottom (2 )  <  0 
ductbottom (2 )  =0; 

end 
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%  This  will  assign  a  true/false  value  if 
%  the  lowest  duct  is  "attached"  to  or 
%  "elevated"  from  the  surface 

if  [ducttop(l)  >  0]  &  [ductbottom ( 1 )  >  hght(l)] 
sf cduct=0 ; 

elseif  ([ducttop(l)  >  0]  &  [ductbottom  ( 1 )  ==  hght(l)])  (*■.«* 
([ducttop(l)  >  0]  &  [ductbottom ( 1 )  ==  0]) 
sf cduct=l ; 

else 

sf cduct=NaN; 

end 

%  Assigns  needed  data  to  the  Plotdata  matrix 
Plotdata (bigcounter, 1) =station; 

Plotdata (bigcounter, 2) =date; 

Plotdata (bigcounter, 3) =time; 

Plotdata (bigcounter, 4) =Ttop; 

Plotdata (bigcounter, 5) =Dewtop; 

Plotdata (bigcounter, 6) =ducttop (1) ; 

Plotdata (bigcounter, 7 ) =ductmiddle ( 1 )  ; 

Plotdata (bigcounter, 8 ) =ductbottom ( 1 ) ; 

Plotdata (bigcounter, 9) =ducttop (2) ; 

Plotdata (bigcounter, 10) =ductmiddle (2)  ; 

Plotdata (bigcounter, 1 1 ) =ductbottom (2 ) ; 

Plotdata (bigcounter, 12 ) =ducttop ( 3 ) ; 

Plotdata (bigcounter, 13) =ductmiddle (3)  ; 

Plotdata (bigcounter, 14 ) =ductbottom ( 3 ) ; 

Plotdata (bigcounter, 15) =sfcduct; 

Plotdata (bigcounter, 16) =Presstop; 

Plotdata (bigcounter, 17) =Pressmid; 

Plotdata (bigcounter, 1 8 ) =Tmid; 

Plotdata (bigcounter, 19) =Dewmid; 

Plotdata (bigcounter, 20) =Mtop; 

Plotdata (bigcounter, 2 1 ) =Mmid; 

Plotdata (bigcounter, 22 ) =pres ( 1 ) ; 

Plotdata (bigcounter, 23) =drct (1) ; 

Plotdata (bigcounter, 24 ) =sknt ( 1 )  ; 

Plotdata (bigcounter, 25 ) =ducttop ( 1 ) -ductbottom ( 1 ) ; 

Plotdata (bigcounter, 26) =ducttop (2 ) -ductbottom (2 ) ; 

Plotdata (bigcounter, 27 ) =ducttop ( 3 ) -ductbottom ( 3 ) ; 

Plotdata (bigcounter, 2  8 ) =Mtop2 ; 

Plotdata (bigcounter, 2  9 ) =Mmid2 ; 

Plotdata (bigcounter, 30) =Mtop3; 

Plotdata (bigcounter, 31) =Mmid3; 

Plotdata (bigcounter, 32) =Mmid-Mtop; 

Plotdata (bigcounter, 33) =Mmid2-Mtop2 ; 

Plotdata (bigcounter, 34 ) =Mmid3-Mtop3 ; 

Plotdata (bigcounter, 35 ) = (Mmid-Mtop) / (duct top ( 1 ) -ductmiddle ( 1 ) ) ; 
Plotdata (bigcounter, 36) = (Mmid2-Mtop2 ) / (duct top (2 ) -ductmiddle (2 ) ) 
Plotdata (bigcounter, 31) = (Mmid3-Mtop3 ) / (duct top ( 3 ) -ductmiddle ( 3 ) ) 
Plotdata (bigcounter, 38 ) =M ( 1 )  ; 

if  Plotdata (bigcounter, 15) ==1  %finds  M  deficit 
Plotdata (bigcounter, 39)  =  M(l)-Mtop; 

else 

Plotdata (bigcounter, 39)  =  NaN; 

end 
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%Check  for  the  effective  "end-of-f lie" 
if  line  ==  -1 
break 

end 

end 

fclose(fid  in);  %  -  end  of  sounding  "decoder"  portion 


%  prepare  the  name  for  the  ".mat" 

%  Logic:  find  the  "/"s  and  s 

name . 

o, 

o 

indexl=find ( filename  ==  '/'); 
if  length ( indexl )  ==  0 
f irst=l ; 

else 

iend=length ( indexl ) ; 
first=indexl (iend) +1; 

name 

end 

index2=find ( filename  ==  .  '); 

Iast=index2-1 ; 

matf ile=f ilename (l,first:last)  ; 
%clear  indexl  index2  first  last  ie: 

O, 

o 

eval ( [ ' save  '  matfile]  ) 

o, 

o 

disp('".mat"  file  has  been  written 
toe 


lie 

we  can  use  the  first  part  of  the 

%  find  the  "/"  in  the  name 
%  if  no  "/"  in  name 

%  the  last  "/" 

%  first  character  in  matfile 

%  use  characters  before  the  " . " 
ans  filename  jo  iline 

to  your  current  directory.') 
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B.  MATLAB  PROGRAM  TO  CALCULATE  MODIFIED  REFRACTIVITY 


function  [M,N]  =  m  n  profile (pres, tmpc, relh, hght) 

Q, 

O 

%  Purpose:  Calculates  the  vertical  profiles  of  Refractive  Index,  N,  and 
%  the  Modified  Refractive  Index,  M,  for  input  vertical 

%  profiles  of  Pressure,  Temperature  Relative  Humidity  and 


Height . 

o, 

o 

%  Input: 

g, 

o 

o, 

o 


o, 

o 

o, 

o 


%  Output : 
units ) 


pres  =  pressure  (millibars) 
tmpc  =  air  temperature  (Celsius) 
relh  =  relative  humidity  (percent,  %) 
hght  =  height  (meters) 

M  =  Modified  Refractive  Index  (dimensionless  refractivity 


o, 

o 

o, 

o 

o, 

o 

g, 

o 

g, 

o 

o, 

o 

o, 

o 

o, 

o 

g, 

o 

o, 

o 

o, 

o 

o, 

o 


N  =  Refractive  Index  (dimensionless  refractivity  units) 

Local  Variables: 

tmpk  =  air  temperature  (Kelvin) 

e_s  =  saturation  vapor  pressure  (millibars) 

ee  =  vapor  pressure  (millibars) 

References 

Bean  and  Dutton,  1968,  Equation  1.16.  ...  for  N 

Patterson,  et  al .  (1994),  Equation  5,  p.  9.  ...  for  M 

Bolton,  Monthly  Weather  Review,  1980.  ...  for  saturation  vapor 


pres 

%  Huschke,  Glossary  of  Meteorology,  1959,  p.  477.  ...  for  vapor 

pressure 

o, 

o 

%  History 

%  Version  1.0  17  July  2001 

%  Mary  S.  Jordan 

%  Dept,  of  Meteorology,  Naval  Postgraduate  School,  Monterey,  CA 

g,  _ 

o 


%  ...  convert  to  temperature  to  Kelvin 

tmpk  =  tmpc+273.15; 

%  ...  compute  saturation  vapor  pressure,  e_s 

e_s  =  6 . 1 12 *exp ( ( 17 . 67  .*  tmpc) . / (tmpc+243 . 5) ) ; 

%  ...  compute  vapor  pressure,  ee 

ee  =  (relh  ./  100.)  .*  e_s; 

%  ...  compute  Refractive  Index,  N,  and  Modified  Refractive  Index,  M 

o, 

o 

N  =  (77 . 6*pres . /tmpk)  -  ( 5 . 6*ee . /tmpk)  +  ( 3 . 75e5*ee . / ( tmpk . A2 ) ) ; 

M  =  N  +  157*hght . / 1 000 ; 

%  - end  of  function - 
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c. 


MATLAB  POGRAM  TO  CALCULATE  POTENTIAL  TEMPERATURE 
AND  SPECIFIC  HUMIDITY 


function  [thta,q]  =  theta_q (pres, tmpc, relh) 

O, 

o 

%  function  [thta,q]  =  theta_q (pres , tmpc, relh) 

O, 

o 

%  written  by:  Mary  Jordan,  NPS  Meteorology  Dept,  10/9/96 

o, 

o 

%  Purpose:  Calculates  Potential  Temperature  (thta)  and 
%  Specific  Humidity  (q) 

O, 

o 

%  Reference:  Atmospheric  Science  by  Wallace  &  Hobbs,  1977, 
Press 


o, 

o 

o, 

o 

o, 

o 

o, 

o 

o, 

o 

o, 

o 

o, 

o 

o, 

o 

g, 

o 


Input:  pres  =  Pressure  (mb) 

tmpc  =  Air  Temperature  (C) 
relh  =  Relative  Humidity  (%) 

Output:  thta  =  Potential  Temperature  (K) 

q  =  Specific  Humidity  (g/kg) 


Academic 


O, 

o 

g. 

o 


o, 

o 


define  local  constants  for  thermodynamic 
p0=1000. ; R=2  87 . ; cp=1004 . ; L=2 . 5e6; 


equations : 


tmpk  =  tmpc+273 . 155;  %  convert  to  tmpc  to  Kelvin 

O, 

o 

e  s  =  6 . 107 8 *exp ( ( 17 . 2 6939 . *tmpc) . / ( tmpc+237 . 3 ) ) ;  %  saturation  vapor 


pressure 

ee  =  (relh  ./  100.)  .*  e_s; 

g, 
o 

w  =  0.622  *  (ee  ./  (pres  -  ee) ) ; 
q  kg  =  w  ./  (1.0  +  w); 

(kg/kg) 

O, 

o 

q=q_kg*1000; 

O, 

o 

thta  =  tmpk  .*  ( (pO  ./  pres ) . A (R/cp) ) 

(K) 


%  vapor  pressure 

%  mixing  ratio  (kg/kg) 

%  specific  humidity 

specific  humidity  (g/kg) 

%  Potential  Temperature 


END 
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D. 


MATLAB  CODE  TO  FIND  STATISTICS 


x=length (data ( : , 1 ) ) ;  %  number  of  rows 

templ=0 ; %temperory  variables 

temp2=0 ; 

temp3=0 ; 

tmpl=0 ; 

tmp2=0 ; 

tmp3=0 ; 

tempminl=0 ; 

tempmin2=0 ; 

tempmin3=0 ; 

temp=0 ; 

%finds  elevation  and  sf  duct  frequencies 
counter  =  1; 
ductsl=0; 
for  i=l : x 

if  data (counter, 6) >0 
ductsl=ductsl+l ; 

end 

counter=counter+l  ; 

end 

percentlduct= (ductsl/x)  *100; 
table= [ ] ; 

table (1,1) =percentlduct; 
counter  =  1; 
ducts2=0 ; 

for  i=l : x 

if  data (counter, 9) >0 
ducts2=ducts2+l ; 

end 

counter=counter+l  ; 

end 

percent2duct= (ducts2/x)  *100; 
table (1,2) =percent2duct; 
counter  =  1; 
ducts3=0 ; 

for  i=l : x 

if  data (counter, 12 ) >0 
duct s3=ducts 3+1 ; 

end 

counter=counter+l  ; 

end 

percent3duct= (ducts3/x)  *100; 
table (1,3) =percent3duct; 
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%finds  frequency  of  day  ductl 

counter=l ; 

dayductl=0 ; 

day=0 ; 

for  i=l : x 

if  data (counter, 3) ==1200 
day=day+l ; 

end 

if  data (counter, 3) ==1200  &&  data (counter, 6) >0 
dayductl=dayductl+l ; 

end 

counter=counter+l  ; 

end 

percentlductday= (dayductl /day)  *100; 
table (2,1) =percentlductday; 

%finds  frequencu  of  night  ductl 

counter=l ; 

nightductl=0 ; 

night=0 ; 

for  i=l : x 

if  data (counter, 3) ==0 
night=night+l ; 

end 

if  data (counter, 3) ==0  &&  data (counter, 6) >0 
nightductl=nightductl+l ; 

end 

counter=counter+l  ; 

end 

percentlductnight= (nightductl/night) *100; 
table (3,1) =percentlductnight; 

%finds  frequency  of  day  duct2 

counter=l ; 

dayduct2=0 ; 

day=0 ; 

for  i=l : x 

if  data (counter, 3) ==1200 
day=day+l ; 

end 

if  data (counter,  3) ==1200  &&  data (counter, 9) >0 
dayduct2=dayduct2+l ; 

end 

counter=counter+l  ; 

end 

percent2ductday= (dayduct2 /day) *100; 
table (2,2) =percent2ductday; 

%finds  frequency  of  night  duct2 

counter=l ; 

nightduct2=0 ; 

night=0 ; 

for  i=l : x 

if  data (counter, 3) ==0 
night=night+l ; 
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end 

if  data (counter, 3) ==0  &&  data (counter, 9) >0 
nightduct2=nightduct2+l ; 

end 

counter=counter+l ; 

end 

percent2ductnight= (nightduct2/night) *100; 
table (3,2) =percent2ductnight; 

%finds  frequency  of  day  duct3 

counter=l ; 

dayduct3=0 ; 

day=0 ; 

for  i=l : x 

if  data (counter, 3) ==1200 
day=day+l ; 

end 

if  data (counter, 3) ==1200  &&  data (counter, 12 ) >0 
dayduct3=dayduct3+l ; 

end 

counter=counter+l  ; 

end 

percent3ductday= (dayduct3/day)  *100; 
table (2,3) =percent3ductday; 

%finds  frequency  of  night  duct3 

counter=l ; 

night duct 3=0 ; 

night=0 ; 

for  i=l : x 

if  data (counter, 3) ==0 
night=night+l ; 

end 

if  data (counter, 3) ==0  &&  data (counter, 12 ) >0 
nightduct3=nightduct3+l ; 

end 

counter=counter+l  ; 

end 

percent3ductnight= (nightduct3/night) *100; 
table (3,3) =percent3ductnight; 

%find  mean  height  of  first  duct 
counter  =  1; 

j=i; 

for  i=l : x 

if  data (counter, 6) >0 

tempi ( j ) =data (counter, 6) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (4,1) =mean (tempi) ; 
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%find  mean  height  of  second  duct 
counter  =  1; 

j=i; 

for  i=l : x 

if  data (counter, 9) >0 

temp2 ( j ) =data (counter, 9) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table ( 4 , 2 ) =mean ( temp2 )  ; 

%find  mean  height  of  third  duct 
counter  =  1; 

j=i; 

for  i=l : x 

if  data (counter, 12 ) >0 

temp3 ( j ) =data (counter, 12) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (4,3) =mean (temp3) ; 

%  find  medians  of  ducts 
table (6,1) =median (tempi ) ; 
table (6,2) =median (temp2 )  ; 
table (6,3) =median (temp3) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (5,1)  =  median (tempi (find (templ<median (tempi) )) ) 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (7,1)  =  median (tempi (find (templ>median (tempi) )) ) 

%  compute  25th  percentile  (first  quartile)  for  duct2 
table (5, 2)  =  median (temp2 (find (temp2<median (temp2) )) ) 

%  compute  75th  percentile  (third  quartile)  for  duct2 
table (7, 2)  =  median (temp2 (find (temp2>median (temp2) )) ) 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (5, 3)  =  median (temp3 (find (temp3<median (temp3) )) ) 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (7, 3)  =  median (temp3 (find (temp3>median (temp3) )) ) 

%  standard  deviations  of  ducts 

table (8,1) =std (tempi) ; 
table (8,2) =std (temp2 ) ; 
table (8,3) =std (temp3) ; 
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%  find  mean,  median  and  std.  dev.  of  thickness 
%  ductl 
counter  =  1; 

j=i; 

templ=0 ; 
for  i=l : x 

if  data (counter, 25) >0 

tempi ( j ) =data (counter, 25) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (9,1) =mean (tempi ) ; 
table  (11, 1) =median (tempi) ; 
table (13,1) =std (tempi) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (10,1)  =  median (tempi (find (templkmedian (tempi) )) ) 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (12,1)  =  median (tempi (find (templ>median (tempi) )) ) 

%duct2 

counter  =  1; 

j=i; 

temp2=0 ; 
for  i=l : x 

if  data (counter, 26) >0 

temp2 ( j ) =data (counter, 26) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table ( 9 , 2 ) =mean ( temp2 )  ; 
table  (11, 2) =median (temp2) ; 
table (13,2) =std (temp2) ; 

%  compute  25th  percentile  (first  quartile)  for  duct2 
table (10,2)  =  median (temp2 (find (temp2<median (temp2) )) ) 

%  compute  75th  percentile  (third  quartile)  for  duct2 
table (12, 2)  =  median (temp2 (find (temp2>median (temp2) )) ) 

%duct3 

counter  =  1; 

j=i; 

temp3=0 ; 
for  i=l : x 

if  data (counter, 27 ) >0 

temp3 ( j ) =data (counter, 27) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 
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table (9,3) =mean (temp3) ; 
table  (11, 3) =median (temp3) ; 
table (13,3) =std (temp3) ; 

%  compute  25th  percentile  (first  quartile)  for  duct3 
table (10, 3)  =  median (temp3 (find (temp3<median (temp3) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  duct3 
table (12, 3)  =  median (temp3 (find (temp3>median (temp3) ))) ; 

%  frequency  of  surface  duct 
counter=l ; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 7) ==data (counter, 8)  &&  data (counter, 7 ) >0 

temp ( j ) =data (counter, 6) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

a=length (temp) ; 

table (1, 4) = (a/x) *100;  %sf  duct  in  all  day 
table (4,4) =mean (temp) ; 
table  ( 6, 4 ) =median (temp) ; 
table (8,4) =std (temp) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (5, 4)  =  median (temp (find (tempkmedian (temp) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (7, 4)  =  median (temp (find (temp>median (temp) ))) ; 

%  find  mean,  median  and  std.  dev.  of  thickness 
%  surface  duct 
counter  =  1; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 15) ==1  &&  data (counter, 25) >0 
temp ( j ) =data (counter, 25) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (9,4) =mean (temp) ; 
table  (11, 4) =median (temp) ; 
table (13,4) =std (temp) ; 

%  compute  25th  percentile  (first  quartile)  for  surface  duct 
table (10, 4)  =  median (tempi (find (templkmedian (temp) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  surface  duct 
table (12, 4)  =  median (tempi (find (templ>median (temp) ))) ; 
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counter=l ; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter,  3) ==1200  ...  . 

&&  data (counter, 7) ==data (counter, 8)  &&  data (counter, 7 ) >0 

temp ( j ) =data (counter,  7)  ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

a=length (temp) ; 

table (2, 4) = (a/day) *100;  %  sf  duct  in  day 
counter=l ; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 3) ==0 ... . 

&&  data (counter, 7) ==data (counter, 8)  &&  data (counter, 7 ) >0 

temp ( j ) =data (counter,  7)  ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

a=length (temp) ; 

table (3, 4) = (a/night) *100;  %  sf  duct  in  night 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  of  ducts 
%  First  duct 
counter  =  1; 

j=i; 

for  i=l : x 

if  data (counter, 32 ) >0 

tmpl ( j ) =data (counter, 32) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (14, 1) =mean (tmpl) ; 
table  (16, 1) =median (tmpl) ; 
table(18,l)=std (tmpl) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (15,1)  =  median (tmpl (find (tmpl<median (tmpl) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (17,1)  =  median (tmpl (find (tmpl>median (tmpl) ))) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  of  ducts 
%  Second  duct 
counter  =  1; 

j=i; 

for  i=l : x 
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if  data (counter, 33) >0 

tmp2 ( j ) =data (counter, 33) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (14, 2 ) =mean ( tmp2 ) ; 
table  ( 16, 2 ) =median (tmp2 ) ; 
table (18,2) =std (tmp2) ; 

%  compute  25th  percentile  (first  quartile)  for  duct2 
table (15, 2)  =  median (tmp2 (find (tmp2<median (tmp2) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  duct2 
table (17, 2)  =  median (tmp2 (find (tmp2>median (tmp2) ))) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  of  ducts 
%  Third  duct 
counter  =  1; 

j=i; 

temp3=0 ; 
for  i=l : x 

if  data (counter, 34 ) >0 

tmp3 ( j ) =data (counter, 34) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (14, 3) =mean (tmp3) ; 
table  (16,  3) =median (tmp3) ; 
table (18, 3) =std (tmp3) ; 

%  compute  25th  percentile  (first  quartile)  for  duct3 
table (15, 3)  =  median (tmp3 (find (tmp3<median (tmp3) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  duct3 
table (17, 3)  =  median (tmp3 (find (tmp3>median (tmp3) ) ) ) ; 


%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  gradient  of  ducts 
%  First  duct 
counter  =  1; 

j=i; 

tmpl=0 ; 
for  i=l : x 

if  data (counter, 35) >0 

tmpl ( j ) =data (counter, 35) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (19,1) =mean (tmpl) ; 
table  (21, 1) =median (tmpl) ; 
table (23,1) =std (tmpl) ; 
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%  compute  25th  percentile  (first  quartile)  for  ductl 
table (20,1)  =  median (tmpl (find (tmpl<median (tmpl) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (22,1)  =  median (tmpl (find (tmpl>median (tmpl) ))) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  gradient  of  ducts 
%  Second  duct 
counter  =  1; 

j=i; 

tmp2=0 ; 
for  i=l : x 

if  data (counter, 36) >0 

tmp2 ( j ) =data (counter, 36) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (19,2) =mean (tmp2) ; 
table  (21, 2) =median (tmp2) ; 
table (23,2) =std (tmp2) ; 

%  compute  25th  percentile  (first  quartile)  for  duct2 
table (20,2)  =  median (tmp2 (find (tmp2<median (tmp2) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  duct2 
table (22, 2)  =  median (tmp2 (find (tmp2>median (tmp2) ))) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  gradient  of  ducts 
%  Third  duct 
counter  =  1; 

j=i; 

tmp3=0 ; 
for  i=l : x 

if  data (counter, 37 ) >0 

tmp3 ( j ) =data (counter, 37) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (19, 3) =mean (tmp3) ; 
table  (21, 3) =median (tmp3) ; 
table (23,3) =std (tmp3)  ; 

%  compute  25th  percentile  (first  quartile)  for  duct3 
table (20, 3)  =  median (tmp3 (find (tmp3<median (tmp3) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  duct3 
table (22, 3)  =  median (tmp3 (find (tmp3>median (tmp3) ))) ; 
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%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strengt  of  ducts 
%  surface  ducts 
counter  =  1; 

j=i; 

tmpl=0 ; 
for  i=l : x 

if  data (counter, 15) ==1 

if  data (counter, 32 ) >0 

tmpl ( j ) =data (counter, 32) ; 

j=j+i; 

end 

end 

counter=counter+l  ; 


end 

table ( 14 , 4 ) =mean (tmpl ) ; 
table  (16, 4) =median (tmpl) ; 
table (18,4) =std (tmpl) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (15, 4)  =  median (tmpl (find (tmpl<median (tmpl) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (17, 4)  =  median (tmpl (find (tmpl>median (tmpl) ))) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  strength  GRADIENT  of  ducts 
%  surface  ducts 
counter  =  1; 

j=i; 

tmpl=0 ; 
for  i=l : x 

if  data (counter, 15) ==1 

if  data (counter, 35) >0 

tmpl ( j ) =data (counter, 35) ; 

j=j+i; 

end 

end 

counter=counter+l  ; 

end 

table (19,4) =mean (tmpl) ; 
table  (21, 4) =median (tmpl) ; 
table (23,4) =std (tmpl) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (20, 4)  =  median (tmpl (find (tmpl<median (tmpl) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (22, 4)  =  median (tmpl (find (tmpl>median (tmpl) ))) ; 

%  finds  mins  and  maxs  of  first  duct 
counter  =  1; 

j=i; 

tempminl=0 ; 
tempmaxl=0 ; 
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for  i=l : x 

if  data (counter, 8 ) ~=0 

tempminl ( j ) =data (counter, 8) ; 
tempmaxl ( j ) =data (counter,  6)  ; 

j=j+i; 


end 

counter=counter+l  ; 

end 

table (24, 1) =min (tempminl) ; 
table (25, 1) =max (tempmaxl) ; 

%  find  mins  and  maxs  of  second  duct 
counter  =  1; 

j=i; 

tempmin2=0  ; 
tempmax2=0  ; 
for  i=l : x 

if  data (counter, 11) ~=0 

tempmin2 ( j ) =data (counter, 11) ; 
tempmax2 ( j ) =data (counter,  9)  ; 

j=j+i; 


end 

counter=counter+l  ; 

end 

table (24, 2) =min ( tempmin2 ) ; 
table (25, 2) =max ( tempmax2 )  ; 

%  find  mins  and  maxs  of  third  duct 
counter  =  1; 

j=i; 

tempmin3=0 ; 
tempmax3=0 ; 
for  i=l : x 

if  data (counter, 14 ) ~=0 

tempmin3 ( j ) =data (counter, 14) ; 
tempmax3 ( j ) =data (counter,  12)  ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (24, 3 ) =min ( tempmin3 )  ; 
table (25, 3 ) =max ( tempmax3 ) ; 

%  find  mins  and  maxs  of  surface  duct 
counter  =  1; 

j=i; 

tempmin=0  ; 
tempmax=0 ; 
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for  i — 1 : x 

if  data (counter, 15) ==1  &&  data (counter, 8 ) ~=0 

tempmin ( j ) =data (counter, 8) ; 
tempmax ( j ) =data (counter,  6) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (24 , 4 ) =min (tempmin) ; 
table (25, 4) =max (tempmax) ; 

%  finds  mean,  median,  first  and  second  quartiles 
%  std.  deviation  of  Mdeficit  of 
%  surface  ducts 
counter  =  1; 

j=i; 

tempdef icit=0  ; 
for  i=l : x 

if  data (counter, 39) >  0 

tempdeficit ( j ) =data (counter, 39) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (26,4) =mean (tempdeficit)  ; 
table  (28,4) =median (tempdeficit)  ; 
table (30,4) =std (tempdeficit)  ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (27,4)  = 

median (tempdeficit (find (tempdef icit<median (tempdeficit) ) ) ) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (29,4)  = 

median (tempdeficit (find (tempdef icit>median (tempdeficit) ) ) ) ; 

%Finds  mean,  median,  lower  and  higher  quartile  and  std.  deviation  of 
%surface  pressure  if  there  is  a  duct 
counter  =  1; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 6) >  0 

temp ( j ) =data (counter, 22) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (31,1) =mean (temp) ; 
table  (33, 1) =median (temp) ; 
table (35, 1) =std (temp) ; 
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%  compute  25th  percentile  (first  quartile)  for  ductl 
table (32,1)  =  median (temp (find (temp<median (temp) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (34,1)  =  median (temp (find (temp>median (temp) ))) ; 

%Finds  mean,  median,  lower  and  higher  quartile  and  std.  deviation  of 
%surface  pressure  if  there  is  no  duct 
counter  =  1; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 6) ==  0 

temp ( j ) =data (counter, 22) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (31,2) =mean (temp) ; 
table  (33, 2) =median (temp) ; 
table (35,2) =std (temp) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (32, 2)  =  median (temp (find (temp<median (temp) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (34, 2)  =  median (temp (find (temp>median (temp) ))) ; 

%Finds  mean,  median,  lower  and  higher  quartile  and  std.  deviation  of 
%surface  pressure  if  there  is  a  surfaceduct 
counter  =  1; 

j=i; 

temp=0 ; 
for  i=l : x 

if  data (counter, 15) ==1 

temp ( j ) =data (counter, 22) ; 

j=j+i; 

end 

counter=counter+l  ; 

end 

table (31,3) =mean (temp) ; 
table  (33, 3) =median (temp) ; 
table (35, 3) =std (temp) ; 

%  compute  25th  percentile  (first  quartile)  for  ductl 
table (32, 3)  =  median (temp (find (temp<median (temp) ))) ; 

%  compute  75th  percentile  (third  quartile)  for  ductl 
table (34, 3)  =  median (temp (find (temp>median (temp) ))) ; 

%//end  of  file 
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E. 


MATLAB  CODE  TO  FIND  DUCTS  MONTH-BY-MONTH 


%finds  ducts  for  january 
x=length (data ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
monthl= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp=='  01 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

monthl (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  february 

x=length (data ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
month2= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  02 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month2 (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  march 
x=length (data  ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
month3= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  03 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month3 (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  april 
x=length (data ( : ,  1 ) )  ; 
counter  =  1; 
countermonth  =1; 
month4= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp=='  04 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month4 (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 
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end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  may 
x=length (data  (:,!)); 
counter  =  1; 
countermonth  =1; 
month5= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  05 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month5 (countermonth, :)=  data (counter, 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  june 
x=length (data (:,!)); 
counter  =  1; 
countermonth  =1; 
month6= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  06 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month6 (countermonth, :)=  data (counter, 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  july 
x=length (data ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
month7= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp=='  07 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month7 (countermonth, : ) =  data (counter, 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  august 
x=length (data ( : , 1 ) ) ; 
counter  =  1; 
countermonth  =1; 
month8= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  08 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 
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month8 (countermonth,  : ) =  data (counter,  :  )  ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  September 
x=length (data ( : , 1 ) ) ; 
counter  =  1; 
countermonth  =1; 
month 9= [ ] ; 
for  i=l :  x 

temp  =  num2str (data (counter,  2  ))  ; 

tempcpr=  temp== '  09 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

month9 (countermonth,  : ) =  data (counter,  :  )  ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  October 
x=length (data  (:,!)); 
counter  =  1; 
countermonth  =1; 
monthl0= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp== '  10 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

monthlO (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  november 

x=length (data  ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
monthl 1= [ ] ; 
for  i=l : x 

temp  =  num2str (data (counter, 2 )) ; 

tempcpr=  temp=='  11 

if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

monthll (countermonth, : ) =  data (counter, : ) ; 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%finds  ducts  for  december 

x=length (data  ( : , 1 ) )  ; 
counter  =  1; 
countermonth  =1; 
monthl2= [ ] ; 
for  i=l :  x 

temp  =  num2str (data (counter,  2 ))  ; 
tempcpr=  temp== '  12 
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if  tempcpr (2 ) ==1  &&  tempcpr ( 3 ) ==1 

monthl2 (countermonth, : ) =  data (counter, 
countermonth  =  countermonth  +  1; 

end 

counter  =  counter  +  1; 

end 

%//end  of  file 
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